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Het  geluid  van  brekende  golven  weerkaatst  aan  de 
zeebodem.  Door  met  een  verticaal  array  van  hydrofoons  dit 
geluid  te  registreren,  komt  informatie  over  deze  bodem 
beschikbaar.  Kan  hiermee  de  bodemsoort  bepaald  worden? 


t 


Problems  telling 

De  Ned er lands e  marine  heeft  vot^r  he! 
gehruiken  van  sonardata  kennis  van  de 
zeebodem  nodig.  Verse  hi  lie  nde  soonen 
btxiems  hebben  namelijk  heel  andere 
effecten  op  de  voortplanling  van  geluid 
onder  water  en  daarroee  op  de  prestatics  van 
sonarsystemen.  Door  te  meten  hoe  een 
uitgezonden  signaal  na  wee rkaat sing  aan  de 
bodem  aankomt  bij  een  ontvanger,  kunnen 
de  kenmerken  van  die  bodem  bepaald 
worden.  Het  heeft  voor  de  marine  echter 
grote  voordelen  als  hie  r  voor,  in  plaats  van 
een  actief  uitgezonden  signaal.  ruis 
ven>or/aakt  dixir  wind  aan  het  zeeoppervlak 
gebruikt  zou  kunnen  worden.  Dit  heeft 
geleid  tot  de  opdracht  van  de  Nederland se 
marine  aan  TNG  Defensic  en  Veiligheid, 
local ie  Den  Haag,  om  te  onderzoeken  of 
omgevingsrois  (ambient  noise)  gebruikt  kan 
worden  om  snel  de  kenmerken  van  de 
omgeving,  met  name  de  zeebodem,  le 
bepalen.  Dit  onderzoek  is  hiermee  een 


onderdeel  van  Rapid  Environmental 
Assessment,  REA. 

Beschnjving  van  de 
werkzaamhedcn 

Omdat  het  geluid  van  brekende  golven 
voortdurend  verandert  en  niei  op  v*>orhand 
bekend  is,  worth  het  effect  van  weerkaatsing 
ervan  aan  de  zeebodem  ge  meten  door  het 
versehil  te  bepalen  tussen  hot  geluid  dal  van 
bo ven  en  van  onderen  binnenkomt  op  een 
verticaal  *  array'  van  hydrofoons.  Het  verties 
aan  int  ensile  it  door  de  reflectie  aan  de  bodem 
is  afhankelijk  van  de  htjek  van  ontvangst  en 
van  de  frequentie  (de  toonhoogte)  van  de 
ruis4  De  manic r  waarop  de  kenmerken  van 
de  bodem  inwerken  op  dit  patrtxm  van 
bodemverlies  hebben  wij  besehreven  met  een 
model.  De  mate  van  detail  van  dit  model  is 
zo  laag  mogelijk  gekozen,  terwij!  het  patnxrn 
toch  reproduceeibaar  blijft.  De  combinatie 
van  txxJemparameiers  die  het  versehil  tussen 
het  ge  meten  en  het  berekende  patroon  zo 
klein  mogelijk  maakt  gebruiken  wij  als 
schatter  van  de  bodemkenmerken. 

Deze  aanpak  wordt  Matched  Field 
Processing  (MFP)  genoemd,  Ondanks  dat 
de  bodem  met  sleehts  vijf  parameters 
beschreven  wordt,  is  bet  aamal  mogelijke 
combinaties  veel  te  grool  om  *met  de  hand' 
de  beste  te  zoeken.  Daarom  hebben  wij  een 
mlatief  nieuwe,  zeer  effici^nte,  aulomatischc 
zoekmethode  verder  ontwikkeld  en 


onderzocht.  LX  combinalie  van 
bodem  verliesmodel  en  automaiische 
zoekmethode  is  eerst  ioegepast  op 
syntheiische  data  en  daama  op  ge  meten 
data.  Spec  tale  aandacht  hebben  wij  be  steed 
aan  het  slatisdsch  interpteteren  van  de 
gevonden  bmlemkenmerken. 

Resukaten  en  condusies 

Ook  al  zijn  deze  ide£en  ongecompliceerd, 
toch  is  het  op  voorhand  onzeker  of  ze  ook 
werken.  tmemalionaal  is  hier  de  afgelopen 
jaxen  veel  onderzoek  naar  gedaan.  In  plaats 
van  alles  zelf  uil  te  vinden,  hebben  wij  oils 
met  ons  eigen  onderzoek  daar  bij 
aangesloten.  IX*  resultaten  van  ons  onderzoek 
in  het  kader  van  dit  M1T-RK  A -project  zijn 
gepresenteerd  in  anikelen  in  internal ionale 
lijdschriften,  De  daarin  gegeven  literatuur- 
verwijzingen  geven  een  actueel  overzicht 
van  het  <>nderzi>ekslerrein. 

Toepasbaurhciti 

Het  resukaat  van  ons  onderzoek  bevestigl 
dat  -  onder  vtxyrwaarden  -  de  kenmerken 
van  de  zeebodem  inderdaad  met  behulp  van 
windruis  bepaald  kunnen  worden. 

De  vwuitzichten  zijn  zelfe  zo  veclbctovend 
dat  cotnmcrciele  hrma  s  operationele 
systemen  ontwikkelen  bestaande  uit  een 
kort  verticaal  array  van  hydrofoons,  hangend 
aan  een  hoei  die  de  stro ruing  van  de  zee 
volgt.  Verwachl  kan  worden  dat  dergelijke 
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systemen  hinnen  lien  jaar  op  tie  rmrkt 
verse hijnen,  waarmee  de  kennis  van  tie 
zeebodem  aan/.ienlijk  toeneemt  en  daarmee 
de  geh miles mogelijkhe den  van  sonar. 
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1  Introduction 


For  the  assessment  of  sonar  performance  in  shallow  water,  environmental  knowledge  is 
required.  Increased  deployment  of  navy  units  in  shallow  water  areas  therefore  has  made 
Rapid  Environmental  Assessment  (REA)  a  hot  issue.  Combining  this  urgent  need  with 
the  idea  to  make  use  of  ambient  noise  data  is  too  attractive  to  refuse.  Bottom  reflection 
loss  can  be  derived  from  ambient  noise  received  on  a  vertical  array.  By  matching  this 
loss  to  modeled  loss,  a  method  called  Matched  Field  Processing  (MFP),  geoacoustic 
parameters  can  be  estimated.  This  idea  has  been  investigated  in  the  research  project 
MFP-REA.  An  overview  of  the  ambient  noise  research  is  presented  in  this  report. 

The  advantages  of  making  use  of  ambient  noise  to  get  environmental  information  are 
obvious.  No  specialised  sound  source  is  needed  and  during  the  recording  of  the  data 
even  no  ship  is  needed,  allowing  for  covert  operations.  However,  noise  is  an  awkward 
phenomenon,  difficult  to  grasp  theoretically,  A  simple  method  to  handle  it  w'ould  be  a 
considerable  advantage.  Such  a  method  is  found  and  is  the  core  of  this  research  project, 
it  is  based  on  the  idea  of  comparing  the  upward  and  dow:nwrard  received  noise  on  a 
vertical  array. 

No  REA  devices  are  yet  available  for  the  navy's  operations.  Bui  in  depth  research  has 
resulted  in  state  of  the  art  knowledge  of  this  research  area.  And,  what  is  more,  operational 
devices  for  REA  are  maturing  fast. 

This  report  discusses  subsequent  articles  that  present  the  research  that  has  been  performed 
for  the  MFP-REA  project.  A  part  of  these  articles  is  available  in  the  appendix. 

The  first  one  presents  the  idea  by  means  of  the  underlying  physics.  In  the  second  paper 
five  experimental  datasets  are  used  to  elaborate  these  ideas  and  to  discuss  data  quality 
issues.  Environmental  parameters  can  be  derived  by  inversion  methods,  which  require 
efficient  global  search  methods.  A  relatively  newr  method,  called  Differential  Evolution, 
is  show  n  to  outperform  existing  methods  by  an  order  of  magnitude  in  article  number 
three  and  four.  Number  five  discusses  the  limited  information  content  of  reflection  loss 
data  from  ambient  noise.  Paper  six  presents  a  theoretical  analysis  of  the  effects  of 
operational  issues  such  as  tilt  in  a  vertical  array,  while  number  seven  assesses  the 
performance  of  a  prototype  freely  drifting  vertical  line  array  buoy.  An  alternative  method 
to  derive  sediment  thicknesses  directly  from  reflection  loss  data  is  described  in  article 
number  eight,  which  method  is  compared  with  the  inversion  approach  in  the  nineth  paper. 
Number  ten  presents  the  efforts  of  a  commercial  firm  to  develop  a  buoyed  array  for 
rapid  environmental  assessment  by  means  of  ambient  noise. 
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2  Overview  of  publications 


To  avoid  violation  of  copyrights,  not  all  the  articles  summarized  in  the  following 
paragraphs  are  reproduced  in  the  appendix.  The  full  text  of  the  papers  of  paragraph  2 A, 
2,2,  2,3,  2.4  and  2,10  are  in  the  appendices  of  this  report.  For  the  other  papers  we  refer 
to  their  respective  journals ^ 

2.1  Geoacoustic  inversion  of  ambient  noise:  a  simple  method  (1 ) 

C.  H  Harrison,  D.G,  Simons,  Proceedings  of  the  Conference  on  Acoustical 
Oceanography,  Southampton  UK,  9-12  April  2001,  ( Proceedings  of  the  Institute  of 
Acoustics ,  Volume  23,  2001,  pp  91-98). 

Pioneering  work  on  inverting  ambient  noise  to  deduce  geoacoustic  parameters  has  been 
successfully  accomplished.  This  is  done  by  comparing  measured  noise  with  that  derived 
from  a  detailed  noise  model.  This  model  requires  as  input  a  noise  source  distribution  and 
environmental  parameters  such  as  bathymetry*  sound  speed  profile  and  bottom  reflection 
properties.  It  has  been  shown  that,  for  range  and  azimuth  independent  environments  and 
with  uniformly  distributed  noise  sources,  geoacoustic  parameters  can  be  inverted  by 
comparing  measured  and  calculated  noise,  searching  over  a  parameter  space, 

However,  detailed  and  time  varying  input  (such  as  sound  speed  profile)  restricts  the 
applicability  of  inversion  with  such  a  complicated  noise  model  This  article,  therefore, 
proposes  an  alternative  approach,  based  on  the  strong  influence  of  the  seabed  reflections 
on  the  vertical  directionality  of  ambient  sea  surface  noise. 

Instead  of  matching  a  measured  to  a  modeled  noise  signal,  the  proposed  method  matches 
reflection  loss  data.  Seabed  reflection  loss  can  be  derived  relatively  easy  from  noise 
received  on  a  vertical  array,  by  separating  the  upward  from  the  downward  field  by  means 
of  beamforming.  For  a  range  independent  environment  and  with  uniformly  distributed 
sea  surface  noise,  a  simple  model  suffices  to  reconstruct  the  plane  wave  reflection 
coefficient  of  the  seabed.  By  matching  this  'measured'  and  modeled  bottom  loss 
(i.e,  the  plane  wave  frequency  and  angle  dependent  reflection  coefficient  in  dB  units), 
inversion  of  geoacoustic  parameters  is  possible. 

The  simplicity  of  the  method  to  derive  the  seabed  reflection  loss  (or  bottom  loss)  can  be 
understood  by  considering  distant  point  sources,  from  which  rays  extend  to  the  receivers 
after  multiple  surface  and  bottom  reflections.  For  every  ray  arriving  directly  (i.e.  most 
recently)  from  the  surface,  there  is  a  corresponding  ray  arriving  from  the  seabed  with 
exactely  one  extra  bottom  reflection,  the  same  number  of  surface  reflections  and  almost 
the  same  horizontal  wavenumber.  So  the  intensity  of  the  noise  from  all  distant  point 
sources,  regardless  of  their  strength,  will  appear  wreaker  by  the  bottom  reflection 
coefficient  Rt( 6)  for  downward  angle  -8  than  for  equal  upward  angle  +#.  In  addition 
the  "distant  source*  condition  can  be  relaxed,  to  a  source  distribution  that  is  close  to 
uniform. 


i 


For  most  of  these  papers  a  copy  for  personal  use  can  also  he  acquired  within  a  few  days  from  Lhe  library 
of  TNO  Ifcienee,  Security  and  Safety.  Phone  +31  70  37406  13  or  email  aertvisC^tno.nl 
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The  bottom  reflection  coefficient  Rh(0f)  not  only  depends  on  the  elevation  angle  0,  but 
also  on  the  frequency/.  This  frequency  dependence  not  only  comes  from  bottom 
interaction,  but  also  from  the  beamforming  required  to  measure  the  noise  intensity 
arriving  in  the  and  -0  directions  on  the  vertical  array.  For  low  frequencies  typically 
around  100  Hz,  the  angular  resolution  of  the  array  degrades  because  of  the  limited 
length  of  the  array.  At  frequencies  above  the  design  frequency  of  the  array  aliasing 

A- 

distorts  the  reflection  loss  pattern.  Hence,  a  beam-smeared  version  Rh  of  the  true 
reflection  coefficient  Rh  is  obtained. 

Preliminary  application  of  this  method  on  experimental  data  gives  promising  results. 

For  a  moderately  range  dependent  environment,  in  the  absence  of  nearby  shipping  and 
for  a  nearly  perfectly  vertical  array,  the  reflection  loss  for  subsequent  elevation  angles 
can  be  derived  accurately,  especially  for  higher  frequencies.  The  technique  is  simple 
and  requires  little  prior  knowledge.  Its  performance  needs  to  be  demonstrated  on 
experimental  data  under  varying  circumstances.  Topics  for  further  study  are  detailed 
inversion  schemes  to  derive  geoacoustic  parameters  from  the  bottom  reflection  loss 
pattern  and  counteracting  the  beam  resolution  problem  at  low  frequencies. 

2»2  Geoacoustic  inversion  of  ambient  noise:  a  simple  method  (2) 

C.H ,  Harrison ,  £>.G.  Simons,  Journal  of  the  Acoustical  Society'  of  America, 

Volume  112  (number  4),  October  2002,  pp  1377  -  1389. 

In  principle  it  is  possible  to  derive  seabed  properties  from  the  bottom  reflection  toss  of 
noise  received  on  a  vertical  array.  In  this  paper  the  potential  of  this  method  is 
demonstrated  for  experimental  data  from  five  sites  (with  water  depths  from  80  to 
140  m).  Special  attention  is  given  to  data  quality  issues*  while  sophisticated  techniques 
for  searching  parameter  space  are  not  applied.  The  main  point  is  to  demonstrate  that  a 
solution  is  possible. 

The  data  are  aquired  with  an  experimental  system,  consisting  of  a  moored  64  element 
vertical  line  array  (VLA)  roughly  centered  in  the  water  column.  Only  the  middle 
32  regulary  spaced  hydrophones  with  half-meter  separation  were  used  for  the  noise 
calculations.  The  electronics  provided  a  flat  band  between  about  100  Hz  and  the  design 
frequency  of  1500  Hz.  The  sample  frequency  was  6  kHz.  (These  data  are  also  used  in 
die  previous  article.) 

A  theoretical  analysis  shows  that  the  approach  can  tolerate  arbitrary  arrangements  of 
distant  ships  and  simultaneous  local  or  distant  wind  noise,  but  that  data  containing 
nearby  point  sources  that  do  not  average  out  over  time  should  be  removed. 

Reflection  loss  versus  frequency  and  angle  is  derived.  Some  propagation  models  can 
directly  use  this  reflection  loss,  but  if  frequency- independent  geoacoustic  parameters 
are  required,  different  analysis  procedures  can  be  applied.  These  procedures  are  not 
elaborated,  since  the  main  objective  of  this  article  is  to  focus  on  the  derived  (beam- 
smeared)  bottom  reflection  loss  data. 

To  take  the  effect  of  the  finite  beamwidth  on  the  reflection  loss  into  account  a  full  noise 
model  would  be  required.  For  demonstration  purposes  the  derived  (beam-smeared) 
bottom  reflection  loss  data  are  compared  here  with  reflection  toss  predictions  under 
simplified  assumptions.  With  uniformly  distributed  dipoles  as  noise  sources  and  a  range 
independent  environment  the  directional  noise  can  easily  be  calculated  and  the  array 
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response  can  be  estimated,  knowing  the  beam  pattern.  These  calculations  are  based  on  a 
three- layer  (two  boundary)  model  of  the  environment. 

For  several  experimental  datasets  the  up/down  beamformed  reflection  loss  is  calculated 
and  presented  as  an  angle-frequency  color  plot.  For  comparison,  ‘correct'  geoacoustic 
parameters  for  each  specific  site  are  taken  from  several  reference  sources. 

These  parameters  are  slightly  adjusted  to  get  a  better  agreement  between  the  measured 
and  the  calculated  beamformed  reflection  loss  patterns.  The  calculated  patterns  are 
remarkably  similar  to  the  experimental  plots,  demonstrating  the  potential  of  the 
proposed  method. 

For  each  experimental  data  set  the  specific  circumstances  at  the  site  had  to  be  taken  into 
account  for  generating  the  calculated  reflection  loss  picture.  In  one  case  there  was 
hardly  any  wind  during  the  experiment.  In  this  case  anomalous  patterns  appear  in  the 
reflection  loss  picture,  resulting  from  sidelobes  of  the  beamforming  to  the  remaining 
horizontal  noise,  while  there  is  hardly  any  noise  coming  from  up  and  down  directions. 
However,  if  this  happens,  it  is  very  obvious. 

Several  arguments  indicate  that  the  reflection  loss  derived  with  this  method  comes  from 
the  vicinity  of  the  first  bottom  bounce,  which  is  a  function  of  angle.  The  consequence  of 
this  is  that  a  drifting  array  could  use  this  method  to  survey  spatial  bottom  variations, 

The  article  ends  with  an  overview  of  analytical  and  experimental  issues.  Analytical 
matters  are  the  effects  of  powerful  directional  sources,  the  consequences  of  sidelobes 
during  beamforming  and  the  inherent  inversion  ambiguity.  Nonacoustic  noise  and  array 
tilt  are  experimental  issues.  In  spite  of  all  conditions  that  have  to  be  fullfilled,  die  method 
is  still  very  promising  and  is  now  ripe  for  an  automated  search  for  geo-acoustic  parameters. 

2.3  Improved  performance  of  global  optimisation  methods  for  inversion  problems  in 

underwater  acoustics 

Camiel  van  Moll ,  Dick  G.  Simons,  Proceedings  of  the  7th  European  Conference  on 
Underwater  Acoustics,  Delft,  5-H  July  2004. 

Sometimes  knowing  reflection  loss  patterns  does  not  suffice,  but  geoacoustic  parameters 
are  required.  Matched  field  inversion  problems  in  underwater  acoustics  relate  to  searching 
for  that  combination  of  environmental  parameters  that  maximises  similarity  between  a 
modelled  and  a  measured  pressure  field.  For  instance,  inversion  of  ambient  noise  looks 
for  the  best  combination  of  sediment  sound  speed,  sediment  layer  thickness,  bottom 
attenuation  and  subbottom  sound  speed,  that  makes  the  calculated  reflection  loss  similar 
to  the  measured  loss. 

Geoacoustic  models  mostly  are  strongly  nonlinear,  meaning  that  small  parameter  changes 
can  cause  big  changes  of  the  calculated  model  outcome.  As  a  consequence,  the  best 
match  between  model  and  measurements  is  very  unpredictable  and  global  search  methods 
must  be  applied.  In  essence  these  methods  are  random  search  methods:  a  combination 
of  parameters  is  chosen  at  random  and  the  distance  of  the  resulting  model  outcome  to 
the  data  is  measured.  The  best  combination  is  the  one  that  minimises  a  distance  function 
(also  called  energy  function).  To  reduce  the  number  of  combinations  that  must  be 
considered,  a  global  search  method  applies  certain  rules  to  generate  a  new  random 
combination  of  parameters. 
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Simulated  annealing  and  genetic  algorithms  are  well  known  and  frequently  applied 
global  search  methods.  Differential  evolution  (DE)  is  a  relatively  new  method, 
described  in  detail  here,  for  which  a  better  performance  is  claimed.  DE  resembles 
genetic  algorithms  in  the  sense  that  both  improve  a  population  of  vectors  during 
successive  generations.  In  this  article  the  claim  is  investigated  by  means  of  a  complicated 
distance  function  with  many  local  minima  that  the  search  algorithm  has  to  escape  from. 

As  performance  criteria  for  optimisation  methods,  efficiency  and  robustness  are 
considered.  Since  the  forward  model  calculations  take  most  of  the  computation  time 
during  optimisation,  efficiency  is  measured  by  the  required  number  of  these  calculations, 
typically  a  few  tens  of  thousands.  The  stochastic  nature  of  these  methods  causes  that  even 
a  global  search  run  can  get  stuck  in  a  local  minimum.  Robustness  is  defined  as  the 
probability  of  an  optimisation  run  to  end  in  the  global  minimum.  It  is  shown  that  with 
comparably  high  robustness  (a  probability  of  0.99)  DE  requires  only  2400  function 
evaluations  instead  of  the  25600  of  a  genetic  algorithm,  an  order  of  magnitude 
improved  efficiency. 

It  is  well  known  that  the  performance  of  optimisation  methods  depends  on  the  specific 
energy  function  applied.  Therefore  it  is  not  sure  that  DE  is  the  best  search  method  for 
geoacoustic  inversion.  However,  since  the  energy  function  used  has  features  similar  to 
those  known  to  be  present  in  geoacoustic  inversion  -  various  local  minima,  parameter 
coupling  and  a  varying  sensitivity  to  the  parameters,  -  the  comparison  of  the  methods 
suggests  that  DE  is  a  very  good  method  for  geoacoustic  inversion.  Therefore  DE  will  be 
used  as  the  method  for  automated  search  of  geoacoustic  parameters  by  means  of  bottom 
reflection  loss  data  derived  from  ambient  noise. 

Matched  field  inversion  is  an  application  of  inverse  theory.  Addressing  non-uniqueness 
and  stability  (a  small  perturbation  to  the  data  can  lead  to  a  large  perturbation  to  the 
recovered  model)  are  major  issues  in  inverse  theory.  Because  inverse  theory  estimates 
the  best  fitting  model,  the  uncertainty  of  this  estimate  should  be  specified.  These  topics 
are  not  investigated  in  depth  in  this  article. 

2.4  Inversion  of  shallow  water  ambient  noise  data  by  means  of  differential  evolution  as 

a  global  search  method 

Dick  C.  Simons,  Camiei  van  Mali  Chris  H  Harrison,  Proceedings  of  the  2nd  Workshop 
on  Experimental  Acoustic  Inversion  Methods  for  Assessment  of  the  Shallow  Water 
Environment,  Island  of  Ischia,  28-30  June  2004. 

Deriving  geoacoustic  parameters  from  experimental  data  is  far  more  difficult  than  from 
synthetic  data.  With  experimental  data  a  large  part  of  control  over  the  situation  is 
sacrificed,  especially  the  simplified  description  given  by  the  bottom  model  leads  to 
uncontrolled  differences  with  the  real  world  situation.  However,  using  'ground  truth" 
geoacoustic  parameter  values  derived  from  other  sources  (as  in  article  2),  still  gives  some 
control  over  the  situation.  Sacrificing  this  prior  knowledge  leads  to  entering  a  huge 
parameter  space  -  six  dimensional  in  case  the  environment  is  modeled  with  six  parameters 
*  which  is  unsuited  for  investigating  all  parameter  combinations.  Automated  global 
search  methods  are  required.  The  main  disadvantage  of  these  methods  is  the  introduction 
of  ambiguous  solutions;  combinations  of  parameters  that  result  in  calculated  bottom  loss 
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that  resembles  measured  bottom  loss  very  well.  From  inverse  theory  it  is  known  that 
these  ambiguities  are  inevitable.  Differential  evolution  is  used  as  global  search  model. 

This  article  gives  a  detailed  description  of  the  bottom  loss  models  applied  and  of 
differential  evolution.  Experimental  data  from  the  2000  Mapex  experiment  in  the 
Mediterrenean  are  used.  These  data  are  inverted  with  a  single  sediment  model,  i.e. 
a  single  homogeneous  sediment  layer,  overlying  a  uniform  subbottom.  Of  57  runs, 

52  gave  consistent  results,  allowing  for  an  impression  of  the  uncertainty  of  the  parameter 
estimates.  However,  the  bottom  loss  patterns  derived  from  these  parameter  settings  are 
much  more  in  agreement  with  the  measured  bottom  loss  than  that  from  the  outcome  of 
the  hand-searching  technique  of  article  2.  The  improved  performance  of  differential 
evolution  is  demonstrated. 

Since  the  experimental  bottom  loss  pattern  suggests  the  presence  of  a  second  sediment 
layer,  the  model  is  extended  to  two  sediment  layers  overlying  a  subbottom  half-space. 
Inversion  of  synthetic  data  made  with  this  model  was  succesfully  inverted,  although 
also  a  second  ambiguous  solution  was  found.  Inverting  the  experimental  data  with  a  two 
layer  forward  model  is  done  200  times,  from  which  1 14  succesful  (low  energy)  runs  are 
selected.  These  runs  only  show  small  differences  between  the  first  and  second  sediment 
layer,  both  in  agreement  with  the  single  layer  inversion  results.  This  indicates  that  a 
single  layer  model  suffices  in  this  case. 

2.5  Geoacoustic  inversion  of  ambient  noise  reflection  loss  data 

J.  Dettmer,  S.E,  Dosso,  M.K.  Prior,  C.H.  Harrison,  N.R.  Chapman,  Proceedings  of  the 
7th  European  Conference  on  Underwater  Acoustics,  Delft,  5-N  July  2004 

The  information  content  of  reflection  loss  data  from  ambient  noise  is  limited. 

Therefore  a  full  characterisation  of  the  environment  can  not  result  from  these  data. 

With  realistic  synthetic  reflection  loss  data  -  with  similar  characteristics  as  the 
experimental  data  of  article  2  -  the  extent  to  which  layer  properties  can  be  resolved  is 
investigated.  The  minimum  thickness  of  the  top  sediment  layer  and  the  maximum  depth 
penetration  are  estimated. 

To  derive  geoacoustic  parameters,  inversion  is  applied  to  synthetic  ambient  noise  data. 
The  parameter  vector  that  makes  the  calculated  reflection  loss  most  similar  to  the 
synthetic  loss  gives  a  good  estimate  of  the  environment.  However,  an  inversion  problem 
can  not  be  regarded  solved  as  long  as  the  parameter  uncertainties  are  not  determined. 

In  Bayesian  theory,  the  posterior  probability  distribution  (PPD)  is  considered  the 
complete  solution  of  the  inverse  problem.  This  PPD  gives  credibility  intervals  for  each 
of  the  estimated  parameters  and  their  correlation.  In  this  way.  confidence  in  the  inversion 
solution  can  be  expressed. 

The  forward  model  used  during  inversion,  calculates  the  plane  wave  reflection  coefficient 
R  as  a  function  of  angle  and  frequency  for  an  arbitrary  set  of  fluid,  homogeneous, 
horizontal  layers.  This  reflection  coefficient  is  applied  to  ambient  noise  due  to  breaking 
surface  waves  to  model  the  upward  and  the  downward  field  arriving  at  the  vertical 
array,  The  effect  of  beamforming  -  required  to  separate  both  fields  -  is  included,  which 
smears  out  the  reflection  coefficient.  Finally  Gaussian  distributed  random  noise  is  added 
to  the  synthetic  data. 
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Each  sediment  layer  is  represented  by  4  parameters:  layer  thickness,  sound  velocity, 
density  and  attenuation.  The  minimum  layer  thickness  is  resolved  by  using  a  single  thin 
layer  over  a  half-space.  Several  inversions  with  different  layer  thicknesses  were 
performed.  For  each  thickness  the  standard  deviation  of  the  upper  layer  parameters 
about  their  true  values  is  determined.  For  synthetic  data  with  a  200-1500  Hz  frequency 
range,  optimal  resolution  is  reached  around  0.7  m.  For  thinner  layers,  the  method  cannot 
resolve  the  sound  speed  of  the  layer. 

The  maximum  penetration  depth  is  derived  by  inverting  two  homogeneous  layers  above 
a  half- space.  The  second  layer  is  1  m  thick  (exceeding  the  0.7  m  value)  and  its  parameters 
can  therefore  be  derived  by  inversion.  The  effect  of  increasing  the  thickness  of  the  upper 
layer  hi  on  the  second  layer  parameter  uncertainty  is  investigated.  At  hi  =  5m  the  density 
of  the  second  layer  is  no  longer  resolved  for  the  synthetic  data  under  consideration. 

It  should  be  kept  in  mind  that  frequency  range  and  depth  of  penetration  are  determined 
by  the  array  length  and  hydrophone  density.  However,  this  investigation  has  shown  that 
it  is  possible  to  invert  reflection  loss  data  from  ambient  noise,  where  the  depth  penetration 
is  in  the  order  of  a  few  metres, 

2.6  Geoacoustic  inversion  of  ambient  noise:  a  modelling  study 

M.K.  Donnelly,  G.A ,  Matt,  UDT  Europe  2003 

SEA,  Systems  Engineering  &  Assessment  Ltd,  is  undertaking  research  with  a  dedicated 
ambient  noise  geoacoustic  inversion  system,  based  on  a  drifting  VLA  buoy.  For  this 
system  a  theoretical  study  aims  at  understanding  key  issues  affecting  performance. 

The  system  consists  of  a  32  hydrophone  VLA  suspended  underneath  a  buoy.  The  design 
frequency  of  the  array  is  4  kHz,  ke.  a  hydrophone  spacing  of  0.1875  m.  The  sampling 
frequency  is  12.288  kHz,  In  this  way  most  of  the  frequencies  of  interest  are  covered. 
(Notice  that  the  frequency  range  of  this  prototype  system  has  shifted  to  frequencies 
higher  than  those  of  the  experimental  system  of  article  number  2.  ),  The  environment 
taken  into  consideration  is  a  range-independent  100  m  deep  isovelocity  water  layer 
overlying  a  two  layer  seabed. 

For  several  situations  synthetic  data  are  generated,  by  successively  modelling  the 
vertically  directed  noise  field  for  clearly  defined  environments  (under  simple  conditions), 
calculating  the  response  of  the  VLA  (without  array  shading)  and  deriving  the  seafloor 
reflection  coefficient.  On  such  data  a  reflection  coefficient  analysis  algorithm  is  applied, 
which  resolves  the  first  sediment  layer  thickness  and  sound  speed,  based  on  interference 
fringes  which  occur  at  high  elevation  angles.  These  parameter  values  are  used  as  initial 
estimates  for  the  inversion  algorithm. 

Theory  states  that  the  down-to-up  in-beam  noise  power  ratio  resolves  the  seafloor 
reflection  coefficient  to  within  the  angular  resolution  limits  of  the  heamformer. 

This  holds  perfect  for  range-independent  environments  w'ith  a  uniform  sea  surface  noise 
source  distribution.  Modelled  results  indicate  that  the  method  still  works  with  moderate 
range  dependent  environments  and  point  sources  in  the  far  field.  This  sets  the  stage  for 
three  performance  issues  that  need  closer  inspection.  What  is  the  effect  of  VLA  tilt? 
What  is  the  effect  of  an  abrupt  seabed  change  at  a  given  range  on  the  performance  of 
the  system?  And  what  effect  has  a  nearby  point  source  (i,e,  ship)? 
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Analysis  of  the  consequences  of  array  tilt  shows  that  moderate  tilt  angles  of  up  to  5°  do 
not  introduce  significant  errors.  The  rule  of  thumb  is  that  errors  introduced  by  tilt  angles 
less  than  the  beamwidth  are  tolerable. 

With  an  abrupt  change  in  the  thickness  of  the  first  sediment  layer  from  1  in  to  3  m  at  a 
range  of  50  m  the  reflection  coefficient  analysis  algorithm  still  gives  acceptable 
estimates  for  sediment  thickness  and  sound  speed  of  the  first  layer  Since  this  algorithm 
makes  use  of  the  fringe  pattern  at  steep  angles,  it  inherently  provides  high  spatial 
resolution  and  resolves  the  seafloor  close  to  the  VLA, 

A  strongly  non-uniform  noise  source  distribution  violates  the  assumptions  on  which  the 
method  is  based.  With  a  ship  at  a  range  up  to  3  km  the  method  breaks  down.  But  for  the 
ship  at  5  km  distance  the  analysis  algorithm  still  gives  acceptable  results  (at  100  m 
water  depth). 

Attention  must  be  given  to  the  fact  that  the  analysis  is  based  on  the  outcomes  of  the 
reflection  coefficient  analysis  algorithm  and  not  on  geo-acoustic  parameter  estimates 
resulting  from  full  inversion.  Analysis  of  the  inversion  results  belongs  to  the  next  phase 
of  work. 

2.7  Measuring  seabed  reflection  loss  from  ambient  noise:  an  experiment  in  the  north¬ 

west  approaches  to  the  UK 

Marcus  Donnelly,  Proceedings  of  the  7th  European  Conference  on  Underwater 
Acoustics,  Delft ,  5-8  July  2004 . 

A  data  collection  trial  was  conducted  west  of  Scotland,  a  region  for  which  several 
reference  data  (a  borehole  and  seismic  data)  are  available.  The  freely  drifting  VLA 
buoy  mentioned  in  article  6  had  the  centre  of  the  array  at  58  m  depth,  approximately 
halfway  the  water  column  of  137  m.  Both  ambient  noise  and  noise  from  the  trials 
vessel  was  recorded. 

The  seabed  reflection  coefficients  derived  from  measured  ambient  noise  are  compared 
w'ith  modeled  reflection  loss  from  the  reference  acoustic  data.  The  uncertainty  margin 
for  these  reference  data  caused  significant  variation  in  the  modelled  reflection 
coefficients.  The  measured  bottom  loss  from  the  ambient  noise  data  agreed  well  with 
one  of  the  modelled  reflection  loss  figures,  but  this  is  no  decisive  confirmation  that  the 
so  called  1  measured  bottom  loss*  is  reliable. 

As  an  extra  check  of  this  4 ambient  noise  bottom  loss*,  it  is  used  to  calculate  the  modelled 
transmission  loss  over  increasing  ranges  of  the  sound  coming  from  the  trial  vessel. 
Because  this  sound  was  recorded  at  increasing  ranges,  the  measured  transmission  loss  is 
known.  Both  are  broadly  in  good  agreement.  Also  the  modelled  transmission  losses  are 
calculated  tor  three  seabed  parameter  combinations  that  are  all  within  the  reference 
uncertainty  bounds.  The  transmission  loss  derived  from  the  ambient  noise  bottom  loss 
agrees  very  well  with  the  modeled  transmission  loss  for  the  *High  reflectivity*  reference 
case. 

From  this  investigation  it  was  concluded  that  bottom  reflection  loss  derived  from  ambient 
noise  data,  collected  with  a  relatively  short  (6  m)  VLA  suspended  from  a  freely  drifting 
buoy,  can  provide  useful  and  reliable  information  about  the  seabed. 
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2,8  Sub-bottom  profiling  using  ocean  ambient  noise 

Chris  Harrison.  Journal  of  the  Acoustical  Society  of  America,  Volume  / 15 
(number  4),  April  2004t  ppI505-15l5 , 

Sound  reflecting  off  a  boundary  experiences  a  reduction  in  amplitude  as  well  as  a  phase 
shift.  Both  effects  are  modeled  in  the  reflection  coefficient.  The  bottom  reflection 
coefficient  inferred  from  directional  ambient  noise  measured  with  a  VLA*  is  derived  as 
a  ratio  of  downward  to  upward  received  noise  intensity.  As  a  ratio  of  intensities 
(incoming  sound  energy  per  square  meter  per  second)  this  reflection  coefficient  only 
contains  amplitude  information.  What  is  normally  called  the  reflection  coefficient,  is  in 
fact  only  its  amplitude. 

The  missing  phase  information  of  the  reflection  coefficient  can  be  restored  from  its 
amplitude,  under  the  condition  that  for  a  layered  seabottom  the  first  one  of  several 
arrivals  of  noise  reflected  off  multiple  boundaries  is  also  the  strongest  one. 

The  density  difference  between  water  and  sediment  normally  takes  care  of  this  condition. 
The  phase  of  the  reflection  coefficient  can  then  be  reconstructed*  using  a  technique 
called  spectral  factorization. 

Consider  sound  reflecting  off  the  sea  bottom  under  a  single  angle.  In  underwater 
acoustics  its  reflection  loss  is  usually  treated  as  a  frequency  domain  quantity. 

The  physical  processes  of  sound  reflecting  at  multiple  boundaries  can  also  be  described 
in  the  time  domain,  giving  explicit  arrivals  from  seabottom  layers.  Both  descriptions  are 
related  one-to-one  by  a  Fourier  transformation.  The  Fourier  transform  of  the  (frequency 
domain)  reflection  coefficient  (consisting  of  amplitude  and  phase)  gives  the  (time 
domain)  impulse  response  of  the  seabed.  This  impulse  response  shows  time  arrivals  of 
the  noise  reflected  off  subsequent  boundaries  in  the  seabed.  Multiplying  these  time 
arrivals  by  the  sound  speed  in  water  (1500m/s)*  instead  of  the  true  but  unknow  n  layer 
sound  speeds*  gives  approximate  depths  of  the  subsequent  layers. 

The  wray  noise  reflects  off  a  boundary  depends  on  its  angle  of  incidence.  Increasing  angles 
show  increasing  depths  of  the  sea  bottom  layers,  which  can  not  be  true.  It  appears  that 
only  the  higher  grazing  angles*  say  70  -  90  degrees,  give  stable  and  reliable  depths. 

To  derive  a  clear  angle-depth  picture*  the  amplitude  of  the  reflection  coefficient  has  to 
be  supplemented  by  its  phase*  since  w  ithout  the  phase  information  the  order  of  the 
subsequent  boundaries  would  be  ambiguous. 

Apart  from  an  extensive  explanation,  this  article  describes  the  application  of  these 
methods  on  experimental  plane  wave  reflection  coefficients  derived  from  ambient  noise 
arriving  on  a  VLA.  The  results  look  promising.  With  a  drifting  VLA  buoy  it  was 
possible  to  construct  a  bottom  profile  over  5  miles  that  agrees  to  some  extent  with  a 
sub-bottom  profile  obtained  writh  standard  seismic  equipment. 
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2.9  Measurement  of  bottom  Layering  and  reflection  using  ambient  noise 

Chris  H.  Harrison ,  Proceedings  of  the  International  Conference  on  Underwater 
Acoustic  Measurements:  Technologies  and  Results,  Heraklion,  Crete . 

28  June  - 1  My  2005 

This  well  written  article  gives  a  more  extensive  explanation  of  spectral  factorization  and 
investigates  with  simulated  data  the  conditions  for  spectral  factorization.  The  method 
performs  well  if  the  minimum  phase  assumption  is  met.  which  requires  a  high  impedance 
contrast  at  the  water  sediment  interface  and  relatively  low  contrast  at  subsequent  layer 
boundaries. 

More  important  is  the  follow' tag  fundamental  observation.  To  discriminate  between 
upward  and  dow  nward  sea  surface  wind  noise,  a  vertical  line  array  is  needed.  By  means 
of  beamforming  the  reflection  coefficient  of  the  sea  bottom  for  subsequent  angles  is 
derived.  The  more  shallow'  grazing  angles  contain  information  about  the  upper  sediment 
layers,  but  information  about  the  deeper  layers  is  contained  in  the  reflection  coefficients 
at  steeper  angles.  However,  every  vertical  array  has  the  narrowest  beams,  and  therefore 
the  best  resolving  power,  in  the  horizontal  direction.  At  sleep  angles  *  the  end  fire  of  the 
VLA  -  the  interference  fringes  become  blurred  because  of  the  wide  beam  width. 

This  limits  the  depth  range  of  the  method.  Reduction  of  the  beam  width  is  possible  by 
increasing  the  length  of  the  array.  However,  operational  requirements  will  restrict  this 
length.  Increasing  the  number  of  receivers  does  not  result  in  narrower  beams. 

With  this  observation  in  mind,  it  is  clear  that  a  method  that  relies  most  on  the  steep 
angle  reflection  coefficients  not  only  uses  a  small  part  of  the  information  available  in 
the  total  bottom  loss  picture,  but  even  the  pan  that  is  hindered  most  by  the  beam 
smearing  effects  at  steep  angles.  In  this  respect  it  is  sensible  to  consider  the  spectral 
factorization  method  as  a  preprocessing  step  for,  for  instance,  an  inversion  method  that 
makes  use  of  the  full  bottom  loss  information  at  all  grazing  angles. 

2*10  Sub-bottom  profiling  with  ambient  noise  measured  on  a  drifting  vertical  array 

Simons  Dick  G.,  Snellen  Mirjam,  Wetc rings  Ardien  and  Harrison  Chris  H.  Proceedings 
of  the  International  Conference  on  Underwater  Acoustic  Measurements:  Technologies 
and  Results,  Heraklion.  Crete ,  28  June  -  I  July  2005. 

During  an  experiment  in  the  shallow  waters  in  the  Mediterranean  {near  the  Ragusa 
Ridge,  south  of  Sicily),  ambient  noise  is  recei  ved  on  a  drifting  vertical  array. 

From  these  experimental  data  a  time  series  of  bottom  reflection  loss  is  derived  (using 
the  up-down  method),  which  is  input  to  two  different  methods  to  infer  the  layering 
structure  of  the  sea  bottom.  Both  methods  are  compared  by  means  of  their  resulting 
bottom  profiles. 

The  first  method  uses  a  model  to  calculate  the  angle  and  frequency  dependent  reflection 
coefficient  for  a  combination  of  geoacoustic  parameters  of  the  seabottom.  For  each 
geoacoustic  parameter  vector  the  difference  between  the  calculated  and  measured 
reflection  coefficient  can  be  quantified.  Inversion  for  the  best  parameter  vector  (the  one 
with  the  smallest  difference)  is  performed  with  Differentia!  Evolution  as  global  search 
method  (cf  article  4). 
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The  second  method  obtains  the  seafloor  layering  in  a  direct  manner.  The  calculations 
consist  of  supplementing  the  reflection  loss  with  phase  information  and  Fourier 
transforming  the  result*  to  derive  the  time  arrivals  of  noise  reflecting  off  subsequent 
bottom  layers.  Multiplying  the  time  arrivals  by  the  sound  speed  in  water  gives  the  layers 
thicknesses  (cf  article  8). 

During  the  experiment  the  vertical  array  drifted  5  miles  in  12  hours.  With  a  footprint  of 
order  the  water  depth  (l  30m)  each  method  gives  a  considerable  number  of  subsequent 
layer  thicknesses.  As  a  result  one  gets  two  bottom  profiles  each  of  which  can  be  assessed 
for  its  consistency.  Comparing  both  profiles  shows  considerable  agreement  which 
gives  confidence  in  both  methods. 

Although  both  methods  use  the  same  experimental  reflection  coefficients  as  input  data* 
the  time  arrival  method  differs  considerably  from  the  inversion  method  in  several  ways. 
It  derives  its  results  directly  and  therefore  needs  far  less  computations  than  inversion 
does.  Furthermore,  inversion  requires  specification  of  the  number  of  layers  beforehand 
and  it  results  in  estimates  of  seabed  parameters,  such  as  layer  thicknesses,  sound  speeds 
and  densities.  For  the  time  arrival  method  the  number  of  layers  together  with  their 
estimated  thicknesses  are  its  outcome,  without  any  other  parameters.  Since  both  methods 
are  based  on  the  same  data  they  can  supplement  each  other,  with  the  time  arrival  method 
as  a  preprocessing  step  for  the  inversion.  The  random  nature  of  the  inversion  method 
provides  statistical  information,  such  that  the  estimate  for  each  of  the  parameters  is 
accompanied  by  a  confidence  interval. 

2.11  Seabed  profiling  with  ambient  noise 

Marcus  Donnelly,  Sea  Technology,  June  2005. 

The  seabed  properties  can  have  a  significant  effect  on  the  detection  range  of  both  passive 
and  active  sonar  in  shallow  water  environments*  For  sonar  performance  modelling* 
detailed  acoustic  characterisation  of  the  seabed  is  not  required;  it  suffices  to  have  the 
correct  acoustic  reflection  coefficient  as  function  of  grazing  angle  and  frequency* 

In  principle  it  is  possible  to  derive  this  bottom  loss  by  measuring  the  vertical  directivity 
of  sea  surface  noise  received  on  a  VLA.  With  such  an  array*  the  variation  of  noise 
intensity  with  vertical  angle  can  be  measured. 

Based  on  this  idea.  Systems  Engineering  &  Assessment  Ltd.  has  the  intention  to  make 
an  operational  system  for  rapid  environmental  assessment*  As  a  prototype,  a  freely 
drifting  VLA  buoy  has  been  demonstrated  on  three  locations  around  the  United 
Kingdom's  coast.  The  environmental  conditions  during  these  experiments  departed 
from  those  assumed  during  the  set  up  of  the  theoretical  approach.  In  spite  of  this*  the 
derived  reflection  loss  patterns  agree  with  patterns  to  be  expected  from  the  knowm 
seabed  types* 

As  soon  as  such  a  system  has  come  to  maturity*  it  is  a  highly  attractive  device.  It  derives 
its  results  quickly  in  terms  of  equipment  deployability  and  data  acquisition*  It  is  able  to 
resolve  geophysical  variability  nearly  real  time  and  has  minimal  impact  on  the 
environment. 
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3  Conclusion 


Using  sea  surface  wind  noise  to  derive  environmental  knowledge  has  big  operational 
advantages  for  navies.  No  sound  source  is  needed  and  no  additional  operational  units 
are  required  during  the  recording  of  the  data  by  the  drifting  VLA  buoy. 

In  the  proposed  method,  the  angle  and  frequency  dependent  bottom  loss  is  derived  from 
the  ratio  of  the  downward  and  upward  received  noise  intensity.  For  a  range  independent 
environment  and  with  uniformly  distributed  sea  surface  noise,  a  simple  model  suffices 
to  reconstruct  this  bottom  loss.  No  complicated  propagation  or  noise  model  is  required 
and  inversion  methods  can  determine  the  unknown  bottom  parameters. 

Both  requirements  of  the  bottom  loss  model  can  be  mitigated.  The  method  tolerates 
limited  variability  of  the  environment  in  range  and  azimuth.  With  a  VLA.  the  footprint 
of  the  method  is  of  the  order  of  the  water  depth.  The  requirement  of  uniformly  distributed 
sources  is  also  relaxed.  The  method  tolerates  unknown  source  distributions  and  noise 
directionality,  it  is  self  compensating.  Even  distant  point  sources  (ships)  are  acceptable. 
Although  tilt  of  the  vertical  array  causes  smearing  over  the  elevation  angles,  the  method 
allows  for  a  moderate  tilt  that  stays  within  bounds  set  by  the  array's  angle  resolution. 

The  performance  of  the  method  depends  on  the  hydrophone  density  and  array  length, 
which  determine  frequency  range  and  depth  of  penetration.  With  half- meter  separation 
of  the  hydrophones,  the  performance  from  1000  to  1500  Hz  is  good  and  the  method  can 
complement  low  frequency  active  inversion  techniques.  The  method  explicitly  shows 
the  critical  angles  (and  thus  sediment  sound  speeds)  and  an  indication  of  the  number  of 
layers.  It  is  tolerant  to  poor  quality  data,  although  unconelated  noise,  such  as  quantization 
noise,  reduces  the  dynamic  range  of  the  reflection  loss.  Nearby  sound  sources,  however, 
lead  to  anomalous  results  and  powerful  directional  sources  can  cause  problems  due  to 
side  lobe  leakage. 

The  method  can  give  a  good  first  assessment  of  the  seabed.  In  those  cases  where  only 
acoustical  properties  of  the  environment  are  needed,  the  method  is  even  more  promising, 
since  it  directly  supplies  accurate  bottom  reflection  loss  data.  The  strength  of  the  method 
makes  commercial  firms  eager  to  develop  the  required  devices,  suggesting  that  within 
several  years  operational  systems  can  become  available. 
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Geoacoustic  Inversion  of  Ambient  Noise:  A  Simple  Method. 
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Abstract 

The  vertical  directionality  of  ambient  noise  is  strongly  influenced  by  seabed  reflections.  Therefore,  potentially, 
geoacoustic  parameters  can  be  inferred  by  inversion  of  the  noise.  In  this  approach t  using  vertical  array 
measurements,  the  reflection  loss  is  found  directly  by  comparing  the  upward  with  the  downward  going  noise* 
Modelling  and  parameter  searching  are  minimised,  and  the  method  does  not  require  a  detailed  knowledge  of  the 
noise  source  distribution. 

1.  Introduction 

It  is  well  known  that  in  shallow  water  the  coherence  and  directionality  of  ambient  noise  depend  on  the  noise 
source  distribution  and  environmental  parameters  such  as  bathymetry,  sound  speed  profile  and  bottom  reflection 
properties.  In  general,  knowing  the  environment,  one  would  require  a  detailed  model  to  predict  noise. 
Nevertheless  pioneering  work  on  inverting  the  noise  to  deduce  gcoaeoustie  parameters  has  been  successfully 
accomplished  [1,2,12]  in  environments  that  are  range-  and  azimuth-independent,  and  where  one  can  rely  on  the 
sources  being  uniformly  distributed.  In  [2]  a  vertically  separated  pair  of  hydrophones  was  used  to  measure  the 
broad  band  noise  coherence,  which  was  then  compared  with  model  calculations  searching  over  geoacoustic 
parameter  space. 

In  the  Mediterranean,  where  shipping  densities  are  high,  the  uniform  source  disfibution  assumption  is 
frequently  violated.  Furthermore  the  source  distribution  may  change  over  a  period  of  hours,  and  there  may  be 
more  dramatic  changes  associated  with  individual  ships  at  short  ranges.  For  this  reason  an  alternative  method  (an 
extension  of  [12])  is  investigated  in  this  paper.  It  is  showm  below  that  it  is  indeed  possible  to  make  deductions 
about  the  seabed  properties  without  knowing  the  detailed  source  distribution.  In  fact  the  reflection  loss  versus 
angle  can  be  found  directly  by  comparing  the  noise  intensity  arriving  from  equal  up  and  dowrn  elevation  angles.  In 
some  applications,  particularly  at  high  frequencies,  a  solution  in  the  form  of  a  reflection  loss  may  be  more 
appropriate  or  useful  than  the  conventional  gcoaeoustie  parameters.  An  additional  benefit  is  that  wrc  do  not  need  a 
detailed  noise  model;  the  penalty  is  that  wc  need  a  vertical  array  to  resolve  angles.  Here  we  investigate  this 
“inversion”  technique  using  some  experimental  ambient  noise  data  from  south  of  Sicily.  Ground  truth  is 
conveniently  provided  by  some  independent  conventional  geoacoustic  inversion  experiments  carried  out  during 
the  same  trial. 

2,  Theory 

The  array  response  to  a  directional  noise  field  and  noise  coherence  can  be  represented  as  a  w'avc  integral  of  a 
Green’s  function  (Equation  9.14  of  [3])  or  as  an  integral  over  ray  arrivals  [4,  5].  Here  wre  choose  the  latter 
approach  in  order  to  allow'  for  spatially  varying  environments  and  source  distributions  as  has  already  been 
developed  in  [6].  Generally,  given  a  noise  directionality  D(<|>,  8)  (a  function  of  elevation  angle  0  and  azimuth  fy) 
and  the  vertical  array's  beam  pattern  B(0,  05)  for  each  steer  angle  0*,  the  array  response  A(0S)  is  given  by 

A(0S }  =  J  | 6)  B{9, 9S )  cos  GdQdt  ( 1 ) 

It  is  shown  in  [6]  that  the  directionality  D(4>,  0)  can  be  expressed  as  two  terms: 

D(&0)  =  Q(d)S(&d)  (2) 

In  the  equivalent  formula  for  coherence  in  [5,6]  the  beam  pattern  is  simply  replaced  by  a  phase  term.  In  equation 
(2),  S  is  a  complicated  function  of  bathymetry,  sound  speed  and  noise  source  distribution  including  multiple 
arrivals,  whereas  Q  is  simply 
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Q(0)  =  Qxp(-asp ),  6  >0 

*^*(^)e*P(-«(«e -*,))»  0<Q 

where  Rb(0h)  is  the  local  (within  one  ray  cycle  of  the  receiver)  bottom  (power)  reflection  coefficient,  and  0t>  is 
related  by  SnelPs  law  to  0  at  the  receiver  The  exponential  terms  represent  respectively  the  attenuations  along  the 
residual  parts  of  the  ray  directly  from  the  surface  (sr)  and  directly  from  the  bottom  (s^-Sp).  If  we  could  access  the 
directionality  itself,  dividing  the  down  D  by  the  up  D  would  give 


D(^-0) 


777^7  =  Rh  (0* )  exp(-a(s,  -  2  sp )) 
y(+v) 


(4) 


Since  this  is  true  for  all  azimuths  <]>  it  is  also  true  for  the  azimuth  integral  ^D(#,8)d# .  Assuming  that  the  effects 
of  absorption  over  one  ray  cycle  are  smalt  we  obtain 


-0)  d<j> 

- RM  (5) 

+#) 

Given  appropriate  angular  resolution  in  (1)  we  can  actually  measure  Rb  as  a  function  of  0  and  hence  0b  since 
the  beam  response  tends  to  the  noise  directionality  D. 

The  reason  for  this  simplicity,  despite  inclusion  of  spatial  variation  of  the  environment  and  noise  source 
distribution,  can  be  seen  as  follows.  Imagine  a  distant  point  noise  source  from  which  a  single  ray  extends,  after 
multiple  surface  and  bottom  reflections,  to  the  receiver.  For  every  such  ray  arriving  directly  (i.e.  most  recently) 
from  the  surface,  there  is  a  corresponding  ray  arriving  from  the  seabed  with  exactly  one  extra  bottom  reflection, 
the  same  number  of  surface  reflections  and  almost  the  same  horizontal  wavenumber.  Provided  the  environment 
changes  only  slowly  over  one  ray  cycle,  the  slight  shifts  in  position  of  reflection  points  have  no  effect.  So  all 
distant  point  sources,  regardless  of  their  strength,  will  appear  weaker  by  Rh(0b)  for  downward  0  than  for  upward 
0.  In  addition  we  can  relax  the  ‘distant  source’  condition  provided  the  source  distribution  is  close  to  uniform  as  in 
the  subsequent  part  of  [6].  Thus  this  approach  can  tolerate  arbitrary  arrangements  of  distant  ships  and 
simultaneous  local  or  distant  wind  noise. 

For  later  reference  we  note  that  the  function  S  in  (2)  reduces  to 


S(fr  8)  =  sin  8t  /[I  -  R{  {8t  )Rh  (8h )  exp (-asc )]  (6) 

with  R,  and  0t  being  surface  loss  and  surface  angle  respectively,  in  the  case  where  the  noise  sources  arc  uniformly 
distributed  dipoles. 


3*  Experiment 

In  May  1999  a  set  of  experiments  (ADVENT99)  was  carried  out  near  the  Adventure  Bank,  south  of  Sicily.  An 
anchored  64-clemcnt  vertical  array  was  used  for  the  noise  measurements.  The  same  array  was  used  in  conjunction 
wfith  sound  sources  both  fixed  and  towed  by  SACLANTCEISPs  ship,  RV  Alliance,  for  conventional  geoacoustic 
inversion  experiments  (unpublished  [7]).  The  62  m  array  was  centred  in  the  80m  water  depth,  and  although  all  64 
hydrophones  were  available  only  the  middle  32  regularly  spaced  ones  (at  half  metre  separation)  were  used  for 
these  noise  calculations.  Noise  was  sampled  at  6  kHz  on  an  automatically  repeated  sequence  of  two  roughly  15 
second  (90112  samples)  windows  every  four  minutes  (at  0.5  and  3.5  min).  The  electronics  provided  a  flat  band 
between  about  100  Hz  and  2  kHz.  Processing  consisted  of  manually  selecting  portions  without  obvious 
contamination  (unpublished  [8])  and  then  beam-forming  them  in  the  frequency  domain  using  MATLAB.  The 
beam-formed  array  response  can  be  determined  conveniently  by  calculating  the  cross-spectral  density  between  all 
pairs  of  channels,  Cjjff),  using  a  ready-made  MATLAB  routine  including  FFT  averaging  [9],  then  multiplying  by 
the  steer  vectors  Wi(0s,f),  Wj(0s,f)  for  frequency  f  and  summing 

A{0S ,/)  =  £2>i  {8t )  Cy  (/)  w  j  {6, )  (7) 

*  j 

An  example  of  the  beam-formed  spectrum  for  the  middle  32  hydrophones  is  shown  in  Figure  1.  One  can 
immediately  see  that  there  is  more  noise  from  above  (top  part  of  figure)  than  from  below  (bottom  of  figure)  and 
that  the  picture  is  dominated  by  the  poor  angular  resolution  of  the  array  at  low  frequencies.  The  design  frequency 
of  this  middle  part  of  the  array  is  1500  Hz  and  there  is  a  clear  grating  lobe  contribution  arcing  from  +90°  at 
1 500Hz  to  +40°  at  2000  Hz.  An  anti-aliasing  filter  cuts  the  response  off  at  about  2  kHz. 
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Figure  1  Beam-formed  ambient  noise  spectrum  measured  on  a  32-elemcnt  vertical  array. 


The  conventional  geoacoustic  inversion  using  low  frequency  sweeps  (150-850  Hz)  [7]  led  to  the  inferred 
optimum  parameters  at  the  ranges  shown  in  Table  1.  3n  the  current  context  we  regard  these  values  as  “ground 
truth"* 


Parameter 

2  km 

5  km 

10  km 

Sediment  speed  cs  (m/s) 

155918 

1568±7 

1 569±32 

Sediment  thickness  h  (m) 

4.310.7 

5.4±0.7 

6,6+6. 7 

Attenuation  a  (dB/'M 

0.4±0.06 

0.5±0.06 

0.610.20 

Density  p  (g/cm3) 

1.53±0.05 

1 .50±0.2 

1.5810.10 

Sub-bottom  speed  cb  (m/s) 

1637117 

1 74215 1 

1702168 

Table  1:  Ground  truth  -  environmental  parameter  estimates  taken  from  [7], 


4.  Simulation 


Because  of  the  problem  of  angle  resolution  it  is  wise  to  estimate  its  effect  on  the  array  response.  Substituting  (2) 
into  (l)  and  dividing  downward  by  upward  array  response  we  obtain  the  ratio  of  beam  responses  tRt  which  is  an 
approximation  to  Rb 


A{-G,) 

AW 


=  »(*,) 


(8) 


We  note  that,  although  S  does  not  truly  cancel  out,  its  influence  on  the  ratio  *R(0S)  is  weak.  Therefore  wc 
assume  that  an  estimate  of  S  under  the  ideal  conditions  of  a  uniform  noise  source  distribution,  as  in  (6),  will 
suffice.  First  we  take  the  ground  truth  data  (for  2  km)  with  an  assumed  sound  speed  of  1510  ms  1  just  above  the 
sediment,  and  we  calculate  the  bottom  reflection  loss  (RL=10log]O(Rb))  versus  angle  using  SAFARI  [10].  Then 
wc  substitute  this  into  (6)  with  Rt=l  and  asc-*0,  and  wc  calculate  the  array  response  ratio  s*R(0fi)  as  distinct  from 
the  directionality  ratio  in  (5). 

Figures  2(a-d)  show  a  set  of  comparisons  between  '31(0)  and  Rb(0)  at  frequencies  of  100,  300.  500,  1000  Hz. 
One  can  immediately  see  the  smoothing  effect  of  the  beam  pattern.  The  interference  nulls  (caused  by  the 
interference  between  the  upper  and  lower  boundaries  of  the  sediment)  arc  covered  up,  and  more  seriously,  at  low 
frequencies  the  beampattern  averages  so  much  that  the  high  loss  part  of  the  RL  curve  (high  angles)  is  biased 
towards  a  low  value.  The  only  solution  is  to  improve  the  angle  resolution  by,  for  instance,  adaptive  beam  forming, 
or  to  attempt  to  deconvolve  the  known  beam  pattern. 
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Figure2  Reflection  loss  calculated  with  SAFARI  using  'ground  truth*  parameters,  Rb  (dashed),  and  ratio  of  array 
responses,  9t  (solid);  at  (a)  100  Hz,  (b)  300  Hz,  (c)  500  Hz,  (d)  1000  Hz. 


The  ‘before*  and  ‘after'  effect  on  the  complete  reflection  loss  vs  angle  and  frequency  can  be  seen  conveniently  by 
first  taking  the  formulae  in  Section  1.6,3  of  [3]  to  evaluate  the  three  layer  problem,  and  then  applying  the  same 
beam  spreading  to  it.  This  is  shown  in  Figures  3(a,b), 

5.  Analysis  of  results 

We  have  several  options  with  the  noise  data.  We  can  collect  together  the  usable  sections  of  noise  with  the  aim 
of  averaging  out  local  changes  in  directionality.  Alternatively  we  can  look  at  15  second  blocks  of  noise  gathered 
at  different  times  of  day  with  the  aim  of  checking  the  sensitivity  of  the  method  to  changes  in  source  distribution. 
In  a  separate  experiment  in  the  vicinity,  MAPEX2G00  (unpublished  [11]),  noise  measurements  made  with  a  256 
element  horizontal  array  detected  about  14  ships  passing  in  4  hours.  It  is  thought  that  there  were  no  ships  close 
enough  to  spoil  the  ‘distant  ship*  approximation  of  Section  3,  although  clearly  a  non-uniform  and  evolving 
distribution  is  to  be  expected.  Invariance  of  the  quantity  1R(0)  in  (8)  would  provide  some  confidence  in  the  theory. 

Figure  4(a-d)  shows  the  quantity  10  logiOM(0)  calculated  from  the  VLA  array  response  according  to  (8)  for 
four  separate  15  second  blocks  each  about  one  hour  apart.  (The  angle  has  been  corrected  from  that  at  the  receiver 
to  that  at  the  seabed  on  the  assumption  that  the  respective  sound  speeds  are  1505  and  1500  ms'1.)  The  array 
response  is  calculated  for  181  steer  angles  from  the  cross-spectral  density  (128-point  FFT,  sampling  at  6000Hz 
with  averaging  over  the  15  seconds  using  the  MATLAB  function  CSD) 
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Figurc3  Reflection  loss  (a)  calculated  with  3-layer  formula  from  [3],  (b)  including  simulated  beam  forming. 


Although  there  arc  slight  differences  the  main  features  appear  not  to  change  with  time.  In  contrast  the  array 
response  for  these  same  times  (Figure  5(a-d))  does  show  significant  variation  particularly  near  horizontal.  This 
suggests,  as  we  hoped,  that  tR(0)  does  not  depend  on  source  directionality  and  is  simply  a  function  of  the  array’s 
position. 

There  are  some  additional  interesting  features  in  Figure  4.  As  expected  from  Figure  2,  at  low  frequencies  10 
logit, tH{0)  is  flatter  than  the  true  10  logJORb(0).  For  frequencies  between  about  200  Hz  and  1500  Hz  (the  design 
frequency  for  this  vertical  array)  the  variation  with  angle  resembles  theory.  There  is  a  slow  rise  in  loss  up  to  the 
critical  angle  followed  by  a  rapid  rise  to  a  plateau  beyond,  as  we  saw  in  Figure  2.  However  the  peak  is  not  at  90°, 
in  fact  Figure  4(a),  in  particular,  shows  a  hint  of  fringe  patterns  following  a  law  of  the  form  fxsin(0)  =  constant. 
Unfortunately  we  might  expect  such  behaviour  from  modal  phenomena,  array  side  lobe  phenomena  or  even 
sediment  layer  thickness  phenomena,  so  we  cannot  use  it  to  discriminate.  Nevertheless  the  resemblance  to  the 
computed  reflection  loss  versus  angle  and  frequency  shown  in  figure  1.24  of  [3]  and  Figure  3  here  is  quite 
striking. 

The  peak  value  of  10  log|0sJl(0)  is  around  10  dB,  whereas  the  plateau  in  Figure  2  is  at  about  l3dB.  Al  these 
steep  angles  the  Rayleigh  reflection  coefficient  for  a  half-space  simply  depends  on  the  density-sound-speed 
product,  so  there  is  limited  scope  for  altering  it.  For  instance,  if  we  set  pc =2  the  bottom  loss  is  9.54dB.  A  more 
likely  explanation  for  the  discrepancy  is  the  addition  of  non-acoustic  noise  to  the  hydrophones.  This  will  have  the 
effect  of  disproportionately  boosting  the  weak  upward-going  signal,  thus  making  the  bottom  appear  to  be  a 
stronger  reflector.  It  is  also  possible,  though  unlikely  in  these  examples,  that  the  upward-going  signal  is 
contaminated  by  the  side  lobes  of  the  often  much  stronger  horizontal  noise. 

In  this  experiment  the  data  were  transmitted  via  a  radio  buoy,  offset  horizontally  from  the  array  via  a  bottom 
cable  at  a  nominal  distance  of  130  m.  A  potential  complication  is  that  splashes  at  the  buoy  would  appear  as 
spatially  biased  noise.  However  the  direct  elevation  angle  to  the  buoy  is  relatively  small  (17°).  Also  in  Figure  5 
there  is  no  obvious  single  feature  at  a  corresponding  angle  for  all  frequencies,  which  suggests  that  the  buoy  is  not 
a  problem. 
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Figure  4  Inferred  Reflection  Loss  versus  angle  and  frequency  at  four  times  (a,b,c,d)  roughly  one  hour  apart. 


Slight  asymmetry  of  the  peaks  in  Figure  5  will  lead  to  an  anomalous  kink  near  zero  grazing  angle  in  reflection 
loss,  as  seen  in  Figure  4  (particularly  (h)  and  (d)).  A  possible  cause  of  this  frequency  dependent  asymmetry  could 
be  that  there  are  several  distant  ships  contibuting  from  different  azimuths  to  a  slightly  tilted  array.  Thus  spectral 
lines  from  various  azimuths  in  the  horizontal  plane  are  interpreted  as  deviating  from  horizontal  in  a  complex 
manner. 

Finally  we  form  the  average  rclection  loss  (Figure  6a)  by  first  taking  the  power  average  of  the  four  array 
response  examples  (Figure  6b). 
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Figure  5  Beam-formed  Array  Response  versus  frequency  at  four  times  (a,b,c,d)  roughly  one  hour  apart. 


6.  Conclusions 

A  simple  technique  for  extracting  gcoacoustic  parameters  from  ambient  noise  measurements  has  been 
proposed  and  tried  out  with  experimental  data.  Unlike  cohercnce-bascd  methods  that  use  a  hydrophone  pair,  this 
technique  docs  not  require  a  noise  model.  It  can  therefore  still  work  when  the  noise  directionality  and  bathymetry 
are  unknown.  The  method  works  best  at  high  frequencies  because  of  the  vertical  array's  limited  angular  resolution 
at  low  frequencies.  Interestingly  this  complements  conventional  active  inversion  techniques  where  sound  speed 
fluctuations  impede  model  matching  at  high  frequencies.  It  is  suggested  that  the  invariance  of  the  inferred 
reflection  loss  (Figure  4)  coupled  with  the  variability  of  the  array  response  (Figure  5)  provides  evidence  that  the 
method  is  sound.  It  is  likely  that  some  of  the  anomalies  found  are  caused  by  poor  quality  data.  For  instance,  the 
maximum  value  of  reflection  loss  (at  high  angles)  is  often  not  as  great  as  expected,  but  this  can  be  explained  by 
contamination  with  non-acoustie  noise.  Since  First  draft  better  quality  data  has  been  collected  at  several  sites 
south  of  Sicily  and  near  Elba. 

Topics  for  further  study  might  include  detailed  inversion  schemes  and  possibly  ways  of  counteracting  the 
beam  resolution  problem. 
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Figure  6  Reflection  Loss  (a)  derived  from  the  power  averaged  array  responses  (b)  seen  in  Figs  4  and  5,  Grey 

scales  unchanged. 
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The  vertical  directionality  of  ambient  noise  is  strongly  influenced  by  seabed  reflections.  Therefore, 
potenti allyT  geoacoustic  parameters  can  be  inferred  by  inversion  of  the  noise.  In  this  approach,  using 
vertical  array  measurements,  the  reflection  loss  is  found  directly  by  comparing  the  upward-  with  the 
down  ward- going  noise.  Theory  suggests  that  this  simple  ratio  is*  in  fact,  the  power  reflection 
coefficient — potentially  a  function  of  angle  and  frequency.  Modeling  and  parameter  searching  are 
minimized,  and  the  method  does  not  require  a  detailed  knowledge  of  the  noise  source  distribution. 
The  approach  can  handle  stratified  environments  and  is  believed  to  tolerate  range  dependence. 
Experimental  data  from  five  sites,  four  in  the  Mediterranean,  one  on  the  New  Jersey  Shelf,  are 
described.  Most  of  the  Mediterranean  sites  had  temporally  varying  noise  directionality,  yet  yielded 
the  same  reflection  properties,  as  one  would  hope.  One  site  was  visited  in  conditions  of  very  low 
surface  noise.  This  paper  concentrates  on  an  experimental  demonstration  of  the  feasibility  of  the 
method  and  data  quality  issues  rather  than  automatic  search  techniques  for  geoacoustic  parameters. 
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h  INTRODUCTION 

It  is  well  known  that  in  shallow  water  the  coherence  and 
directionality  of  ambient  noise  depend  on  the  noise  source 
distribution  and  environmental  parameters  such  as  bathym¬ 
etry,  sound- speed  profile,  and  bottom  reflection  properties.  In 
general,  knowing  the  environment,  one  could  use  a  detailed 
propagation  model  to  predict  noise.  Potentially,  one  could 
therefore  deduce  bottom  properties  from  the  noise  without  a 
specialized  sound  source  or  even  a  ship.  Pioneering  work  on 
inverting  the  noise  to  deduce  geoacoustic  parameters  has 
been  successfully  accomplished  [Buckingham  and  Jones 
(1987);  Carbone  et  aL  (1998);  Aredov  and  Furduev  (1994)] 
in  environments  that  are  range-  and  azimuth  independent, 
and  where  one  can  rely  on  the  sources  being  uniformly  dis¬ 
tributed.  Carbone  used  a  vertically  separated  pair  of  hydro¬ 
phones  to  measure  the  broadband  noise  coherence,  which 
was  then  compared  with  model  calculations  searching  over 
geoacoustic  parameter  space. 

In  the  Mediterranean,  where  shipping  densities  are  high, 
the  uniform  source  distribution  assumption  is  frequently  vio¬ 
lated.  Furthermore,  the  source  distribution  may  change  over 
a  period  of  hours,  and  there  may  be  more  dramatic  changes 
associated  with  individual  ships  at  short  ranges.  For  this  rea¬ 
son  an  alternative  method  [an  extension  of  Aredov  and  Fur¬ 
duev  (1994)]  is  investigated  in  this  paper.  A  preliminary  ac¬ 
count  of  this  work  was  given  by  Harrison  and  Simons 
(2001).  It  is  shown  below  that  it  is  indeed  possible  to  make 
deductions  about  the  seabed  properties  without  knowing  the 
detailed  source  distribution.  In  fact,  the  reflection  loss  versus 
angle  can  be  found  directly  by  comparing  the  noise  intensity 
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arriving  from  equal  up  and  down  elevation  angles,  in  some 
applications,  particularly  at  high  frequencies,  a  solution  in 
the  form  of  a  reflection  loss  may  be  more  appropriate  or 
useful  than  the  conventional  geoacoustic  parameters.  An  ad¬ 
ditional  benefit  is  that  we  do  not  need  a  detailed  noise  model* 
but  on  the  other  hand  we  do  need  a  vertical  array  to  resolve 
angles.  Here,  we  demonstrate  the  potential  of  the  method 
with  experimental  data  from  five  sites  in  the  Mediterranean 
and  one  on  the  New  Jersey  Shelf.  A  future  paper  will  discuss 
techniques  for  automatic  geoacoustic  parameter  searching. 

Ground  truth  is  conveniently  provided  by  independent 
conventional  geoacoustic  inversion  carried  out  either  during 
the  same  experiments  or  on  different  occasions.  There  are 
also  a  number  of  cores  and  geoacoustic  models  in  most  of 
the  areas.  Specific  references  will  be  made  later. 

II,  THEORY 
A.  Noise  theory 

It  is  fairly  obvious  that  the  noise  field  depends  on  the 
reflection  properties  of  the  seabed,  but  not  so  obvious  that 
the  plane  wave  reflection  coefficient  is  so  simply  related  to 
the  noise  directionality.  The  array  response  to  a  directional 
noise  field  can  be  represented  as  an  integral  over  ray  arrivals 
[Chapman  (1988),  Harrison  (1996)]  or  as  a  wave  integral  of 
a  Green's  function  [Eq.  (9,14)  of  Jensen  et  at.  (1994)],  The 
ray  approach  for  a  range -dependent  environment  is  spelled 
out  in  the  Appendix.  Here,  we  avoid  these  cumbersome  deri¬ 
vations  by  choosing  a  simple  energy  flux  argument  to  explain 
the  effect  in  a  more  concise  and  convincing  way. 

Consider  a  block  of  sea  water,  Fig.  1,  bounded  by  the 
sea  surface,  the  seabed,  and  two  vertical  imaginary*  lines 
(planes)  in  an  otherwise  infinite  horizontally  uniform,  verti¬ 
cally  stratified  environment.  Insert  an  infinite  source  layer  at 
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FIG.  I.  Diagram  showing  a  finite  block  from  an  infinitely  wide  range- 
independent  sea.  Although  sound  may  propagate  obliquely  the  horizontal 
flux  balances  out,  so  here  we  only  show  the  vertical  components.  Up-  and 
down-going  components  of  (plane  wave)  flux  ,  SD  emanate  from  a  hori¬ 
zontal  sheet  source  just  below  the  sea  surface.  The  vertical  balance  leaves 
components  V  and  D  in  the  btnly  of  the  water,  which  are  simply  related  to 
Uxe  bottom  reflection  loss  R$  . 


some  small  (finite)  depth  below  the  surface.  Now,  consider 
the  balance  of  the  upward  and  downward  component  of  flux 
in  the  box  (both  above  and  below  this  source  sheet).  Note 
that  because  of  the  horizontal  uniformity  there  is  no  net  hori¬ 
zontal  flux,  so  we  can  ignore  the  horizontal  component.  We 
regard  the  infinite  source  as  emitting  an  incoherent  plane 
wave  (as  a  distant  light  bulb  would!)  at  a  given  angle  which 
subsequently  obeys  Snell's  law  at  the  various  depths.  Let  the 
sheet  source  emit  flux  upwards  and  SD  downwards  (at  the 
given  angle).  Beneath  the  source  we  define  upward  and 
downward  fluxes  (at  the  corresponding  Snell  angle)  as  U  and 
D,  respectively.  Just  below  the  source  the  total  downward 


flux  includes  the  surface  reflected  up  ward- going  flux 

D=Sd+Rs(Su-$-U ),  (I) 

where  Rs  is  the  surface  power  reflection  coefficient. 

By  the  definition  of  the  bottom  power  reflection  coeffi¬ 
cient  Rb  ,  we  have 

U  =  DRh.  (2) 

If  we  define  an  "effective”  source  strength  S  (since  R,  may 
not  be  unity)  as 

S=Sd+RsSu,  (3) 

then  we  find 

D  =  S/(\-RsRh),  (4) 

V  =  RhSfl\-R,Rh).  (5) 


These  are  exactly  the  formulas  obtained  by  Chapman  (1988) 
and  Harrison  (1996)  [also  seen  in  Eqs,  (A6)  and  (A3)  of  the 
Appendix  and  Harrison  (1997a)]  but  with  volume  absorption 
set  to  zero.  Including  absorption  in  the  flux  argument  is 
straightforward.  Dividing  Eq.  (5)  by  Eq.  (4),  we  obtain  the 
up-to-do wn  noise  ratio  which  is  evidently  the  same  as  Eq, 
(2),  This  is  no  surprise  since  it  is  the  definition  of  reflection 
coefficient!  If  we  now  put  a  receiver  in  this  field  we  can 
measure  U  and  />,  and  hence  deduce  Rft ,  but  clearly  the 
original  field  consisted  of  plane  waves  so  we  will  deduce  a 
plane  wave  reflection  coefficient.  Also,  the  ratio  depends  on 
local  bottom  properties,  i.e.,  within  about  one  ray  cycle  for 
the  angle  in  question. 


Note  that  no  matter  where  the  sources  are,  as  long  as 
they  are  far  enough  away  or  large  enough  to  be  considered 
plane  wave  sources,  Eq.  (2)  is  always  true.  Therefore,  the 
underlying  tenet  of  this  paper  still  holds.  The  plane  wave 
condition  reduces  to  being  able  to  ensure  that  f/t  D  at  the 
receiver  are  substantially  the  same  as  at  the  bottom.  Clearly, 
local  point  sources  (that  cannot,  in  some  way,  be  spatially  or 
temporally  averaged  out)  violate  these  assumptions. 

This  is  a  powerful  derivation  for  the  following  reasons: 

(i)  We  have  made  no  assumption  about  rays,  modes,  or 
waves, 

(it)  We  have  made  no  assumption  about  the  dipole  or 
monopole  nature  of  the  source. 

(iii)  We  have  made  no  assumption  about  the  spectrum  of 
the  source. 

(iv)  Refraction  makes  no  difference  to  the  argument  as 
long  as  energy  gets  from  the  surface  to  the  receiver 
via  the  bottom  (see  the  Appendix  for  more  details), 

(v)  A  small  net  horizontal  component  of  flux  does  not 
alter  the  vertical  component.  Therefore,  neither  non- 
uniform  distant  sources  nor  range-dependent  environ¬ 
ments  spoil  the  argument  (see  the  Appendix  for  more 
details). 

Thus,  this  approach  can  tolerate  arbitrary  arrangements  of 
distant  ships  and  simultaneous  local  or  distant  wind  noise. 
However,  nearby  point  sources  (with  slight  asymmetry  in 
eigenray  angles)  that  do  not  average  out  over  time  lead  to 
problems  that  can  nevertheless  be  anticipated  and  isolated  by 
inspection  of  the  array's  beam  response. 

A  dramatic  demonstration  of  the  correctness  of  the 
theory  (at  least  for  range-independent  environments)  is  given 
by  Harrison  and  Baidacci  (2002).  Within  the  confines  of  one 
wave  model,  one  can  construct  both  the  "ground  truth"  re¬ 
flection  loss  and  the  result  of  the  proposed  VLA  processing 
technique.  First,  OASR  [Schmidt  (1999)]  was  used  to  calcu¬ 
late  the  bottom  reflection  loss  for  a  given  environment  (simi¬ 
lar  to  the  Sand  site.  Sec.  IV  A  later);  a  preview  is  shown  in 
Fig.  2(a).  Then,  the  array  response  to  a  sheet  noise  source 
was  found  from  the  correlation  matrix  option  of  oasn 
[Schmidt  (1999)],  and  the  up-to-dowm  noise  ratio  was  ob¬ 
tained  from  it  [Fig.  2(b)],  In  this  example  the  array  spans  the 
full  water  column  and  these  two  plots  are  very  nearly  iden¬ 
tical.  Naturally,  as  the  array  size  is  reduced  the  deduced  re¬ 
flection  loss  gradually  degrades. 


B,  Possible  analysis  procedures 

There  are  several  ways  we  can  process  this  type  of  VLA 
data.  First,  it  is  possible  to  leave  the  output  as  reflection  loss 
versus  frequency  and  angle  since  this  is  exactly  what  is  re¬ 
quired  for  some  propagation  models,  e.g.,  ray  traces.  Further¬ 
more,  it  is  possible  to  make  comparisons  with  other  reflec¬ 
tion  measurement  techniques  such  as  the  "move-out" 
method  where  pulse  travel  times  and  known  geometry  allow 
direct  measurement  [Holland  and  Osier  (2000)],  Second,  it  is 
possible  to  try  to  convert  the  data  to  frequency-independent 
geoacoustic  parameters  by  using  a  search  algorithm  and  a 
plane  wave  reflection  formula  or  model  (as  opposed  to  a 
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noise  model  i.e.,  one  which  needs  to  assume  a  complete 
environment  and  noise  source  distribution — precisely  what 
we  hope  to  avoid  with  this  method).  In  the  present  paper  we 
will  venture  no  further  than  hand  searches,  the  main  point 
being  to  demonstrate  merely  that  a  solution  is  possible.  The 
precise  meaning  of  “hand  search"  will  be  explained  during 
discussion  of  the  examples  in  Sec.  IV — the  point  is  that  there 
is  a  logical  sequence  in  which  one  searches  for  one  param¬ 
eter  at  a  time.  A  third  possibility  is  to  use  amplitudes  and 
fringe  patterns  to  calculate  layer  thicknesses  and  velocities 
directly  by  methods  more  akin  to  those  used  in  geophysical 
prospecting.  By  Fourier  transforming  from  frequency  to  time 
it  is,  in  principle,  possible  to  convert  the  fringes  (being  a 
power  spectrum)  to  the  layers1  autocorrelation  function  (a 
series  of  delta  functions  at  delays  corresponding  to  all  layer 
separations).  This  time-domain  information  is  then  available 
as  a  function  of  angle. 

In  ail  three  methods  we  first  need  to  map  the  up-to-down 
ratio  from  the  angle  measured  (by  beamforming)  at  the  re¬ 


ceiver  6r  to  the  angle  at  the  seabed  Bb .  Thus,  knowing  the 
sound-speed  profile,  by  Snells  law  we  have 

9b  ~  acos(  {cb  /  c  r)  cos(  0T) ) ,  (6) 

and  with  upwrard  refraction,  zero  bottom  grazing  angle  maps 
to  a  finite  nonzero  angle  at  the  array.  Conversely,  w  ith  down¬ 
ward  refraction,  zero  angle  at  the  anay  maps  to  a  nonzero 
angle  at  the  bottom,  and  there  will  be  a  small  range  of  angles 
at  the  bottom  that  is  not  detectable.  It  is  also  necessary  that 
the  path  connect  with  the  sea  surface — otherwise,  there  could 
not  be  any  measured  surface  noise  at  the  given  angle.  If  there 
is  a  sound  speed  higher  than  that  at  the  receiver  depth  any¬ 
where  above  the  receiver  cm,  there  is  a  possibility  of  a 
“noise  notch,”  a  range  of  angles  that  is  surface-noise-free. 
The  edge  of  this  range  obviously  corresponds  to  a  nonzero 
angle  at  the  receiver 

0m  =  acos  (cricm).  (7) 

Beyond  this,  the  bottom  angle  may  be  greater  or  smaller  than 
the  angle  at  the  receiver  according  to  Eq.  (6). 


0  10  30  30  40  SO  60  70  80  00 

Bottom  Angle 
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ft  id  15 

Reflection  Loss  (dB) 


0  10  20  30  4Q  50  60  70  60  90 
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a  to 

Reflection  Loss  (dB) 


15 


FIG,  2,  (a)  Plane  wave  reflection  loss  for  a  hypothetical  throe-layer  envi¬ 
ronment  similar  to  the  Sicily  sand  case  modeled  with  oasr,  and  (h)  the  ratio 
of  upward- to-do wn ward  noise  for  a  0.5-m-spaced  vertical  airay  spanning 
the  water  column  calculated  by  oasn  assuming  a  sheet  source  and  the  same 
environment. 


C.  Simulated  bottom  reflection  and  beamforming 

In  the  subsequent  sections  we  make  comparisons  be¬ 
tween  measurements  and  predictions.  Unless  otherwise 
stated,  the  predictions  are  based  on  a  three-layer  (two- 
boundary')  model  with  the  bottom  layer  supporting  shear 
waves.  In  some  cases  an  extra  intermediate  fluid  layer  has 
been  added.  Theory  is  based  on  the  analytical  solutions  from 
Jensen  et  ai  (1994),  Eqs.  (1.49)  (for  shear),  (1.57)  (for  three* 
layer),  and  (1.65)  (for  extension  to  m  layer).  The  effect  of  the 
finite  beamwidth  could  be  estimated  by  using  a  full-noise 
model  such  as  dasn  [Schmidt  (1999)}  or  canary  [Harrison 
(1997b)],  This  is  taken  up  by  Harrison  and  Baldacci  (2002). 
Instead,  for  demonstration  purposes  we  note  here  that,  in  the 
case  where  the  noise  sources  are  uniformly  distributed  di¬ 
poles  and  the  environment  is  independent  of  range,  the  di¬ 
rectional  noise  can  easily  be  calculated  (see  the  Appendix). 
We  can  therefore  estimate  the  array  response  knowing  the 
beam  pattern.  Dividing  downward  by  upward  array  response 
we  obtain  the  ratio  of  beam  responses  91,  which  is  an  ap¬ 
proximation  to  Rb 

Al  -  0o) 

Thus,  for  each  experimental  plot  we  can  provide  an  unmodi¬ 
fied  layer  model  plot  and  also  a  beam-smudged  version. 


Ill  EXPERIMENTS 

Experiments  have  been  carried  out  at  five  sites  in  the 
Mediterranean  Sea.  Three  sites  are  south  of  Sicily,  near  the 
Ragusa  Ridge  [Fig.  3(a)],  one  near  the  Adventure  Bank  dur¬ 
ing  ADVENT99  in  May  1999,  and  the  other  two  on  the 
Malta  Plateau  during  MAFEX2000bis  in  November  2000, 
The  other  two  are  north  and  east  of  Elba  [Fig,  3(b)]  also 
during  MAPEX2000bi$,  Bottom  types  include  sand,  silt, 
mud,  and  rock,  and  one  of  the  sites  was  visited  twice  first 
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FIG.  3.  (a)  Three  experimental  sites  south  of  Sicily  [Adventure  Bank  (la¬ 
beled  ‘T),  ADVKNT99  and  Ragusa  Ridge  (labeled  “Sa”  and  “ft’*), 
MAPEX2{KX)bis];  (b)  two  sites  to  the  north  ("S”)  and  cast  (“M")  of  the 
Island  of  Elba  (MAPEX2000t>is);  (c)  one  site  on  the  New  Jersey  Shelf 
C*C"),  BOUNDARY 200 1  showing  the  area  extensively  surveyed  by  Goff 
etal  (1999). 


with,  then  without,  significant  wind  noise.  Another  site  was 
visited  on  the  New  Jersey  Shelf  during  BOUNDARY200I 
[Fig.  3(c)]. 

In  all  cases  SACLANTCEN s  moored  64-element  verti¬ 
cal  array  (VLA)  was  used  for  the  noise  measurements 
(Troiano  et  a l .,  1995).  The  same  array  was  used  in  conjunc¬ 
tion  with  sound  sources  both  fixed  and  towed  by  SACLANT¬ 
CEN  Ts  ship,  RV  Alliance,  for  conventional  geoacoustic  in¬ 
version  experiments  (Siderius  et  ui,  2001).  The  62-m  array 
was  roughly  centered  in  the  water  column,  and  although  all 
64  hydrophones  were  available  only  the  middle  32  regularly 
spaced  ones  (i.e.,  total  length  16  m  at  half-meter  separation) 
were  used  for  these  noise  calculations. 

Noise  was  sampled  at  6  kHz  on  an  automatically  re¬ 
peated  sequence.  During  ADVENT99  the  sequence  consisted 
of  two  roughly  15-s  windows  every  4  min  interlaced  with 
active  transmissions.  This  arrangement  led  to  less  than  ideal 
gain  settings  for  the  noise  measurements  and  a  tendency  to 
introduce  quantization  effects.  In  all  the  other  cases  noise 
measurements  were  taken  for  10  s  every  10  s  while  the  RV 
Alliance  was  quiet.  The  electronics  provided  a  flat  band 
between  about  100  Hz  and  2  kHz.  Processing  consisted  of 
manually  selecting  portions  without  obvious  contamination, 
and  then  beamforming  them  in  the  frequency  domain. 

IV,  ANALYSIS  OF  RESULTS 

The  experimental  beamformed  array  response  can  be  de¬ 
termined  conveniently  by  calculating  the  cross-spectral  den¬ 
sity  between  all  pairs  of  channels,  C^(/h  using  a  propri¬ 
etary  routine  including  FFT  averaging,  then  multiplying  by 
the  steer  vectors  wf(0(>,/),  wj( B0  J )  for  frequency/,  and 
summing 

4<0O,/)=WTCW=2  X  <WC;;(/)W;(<U  (9) 

I  j 

Each  array  response  was  calculated  for  181  steer  angles 
from  the  cross -spectral  density  matrix  for  the  middle  32  el¬ 
ements  (17-48)  of  the  64-element  array.  The  bottom  hydro¬ 
phone  (number  I)  was  moored  at  11  m  above  the  seabed  in 
the  Mediterranean  cases  and  14  m  in  the  New  Jersey  case. 
The  cross- spectral  density  used  nonoverlapping  128-point 
FFTs,  sampling  at  6000  Hz,  averaging  over  the  10  s  of  the 
file.  Further  power  averaging  was  carried  out  over  10  to  30 
files.  In  most  cases  a  total  of  about  5  min  of  data  (though 
carefully  selected)  proved  adequate. 

The  Mediterranean  sites  have  fairly  heavy  shipping,  par¬ 
ticularly  those  on  the  Malta  Plateau  and  Adventure  Bank,  In 
a  separate  experiment  in  the  vicinity  (MAPEX2000).  noise 
measurements  made  with  a  256-element  horizontal  array  eas¬ 
ily  detected  about  14  nearby  ships  passing  in  4  h  (Harrison, 
2001).  It  is  thought  that  there  were  no  ships  dose  enough  to 
spoil  the  ‘‘distant  ship”  approximation,  although  clearly  a 
nonuniform  and  evolving  distribution  is  to  be  expected.  At 
the  sites  near  Elba  there  tend  to  be  occasional  ferries  and 
small  fishing  boats  rather  than  container  ships.  In  contrast, 
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TABLE  I.  Important  sound  speeds  for  angle  mapping. 


Depths  (m) 

Sound  speeds  (m/s) 

Min. 

Date 

dd-mm-yyyy 

Site 

Array 

Water  center 

Array  Bottom  Max. 

er  £$  Cjuu 

bottom 
angle  (") 

Wind 

(kt) 

06-05-1 999 

Adv.  Bank 

SO 

39 

1509.4 

1511.2 

1513.5 

5-8 

22*11*2000 

S  Sicily 

139 

96 

1512.0 

1512.0 

1 524.2 

0.00 

20-22 

23*11-2000 

S,  Sicily 

99 

56 

1512.8 

1511.4 

1522.5 

2.47 

5-6 

29*1 1-2000 

N.  Elba 

121 

78 

1520.5 

1512.1 

1520.8 

6.03 

3-5 

30*1 1*2000 

E.  Elba 

125 

82 

1520.4 

1508.8 

1520.6 

7,08 

12-16 

01-12-2000 

N.  Elba 

121 

78 

1519.9 

1513,9 

15199 

5,09 

3-4 

17-05-2001 

N.J.  Shelf 

140 

92 

1499.3 

1499 3 

1500.5 

0,00 

10-14 

there  was  hardly  any  nearby  shipping  at  the  New  Jersey 
Shelf  site.  Wind  speeds  are  shown  in  Table  l  Sound- speed 
profiles  at  the  various  sites  are  shown  in  Fig,  4. 


A.  “Sand1'  site:  S.  Sicily  (near  Ragusa  Ridge),  22 
November  2000 

The  array  response  versus  angle  and  frequency  is  shown 
in  Fig.  5(a).  There  are  various  features  common  to  all  the 
following  cases.  On  the  right  the  intensity  falls  off  rapidly 
because  of  a  2-kHz  antialiasing  filter.  The  general  broaden¬ 
ing  towards  the  left  is  a  manifestation  of  the  poorer  angle 
resolution  at  lower  frequency.  The  double  spikes  near  the 
horizontal  span  a  noise  notch  caused,  as  explained  earlier  by 
the  sound-speed  profile  (see  Fig.  4).  Straight  away  we  can 
see  a  considerable  difference  between  up-  and  down-going 
intensity.  The  "up-to-down”  ratio  is  shown  in  Fig.  5(b). 
Angles  have  been  corrected  from  that  at  the  receiver  to  that 
at  the  seabed  according  to  the  sound  speeds  shown  in  Table  I. 
(Note  that  Table  I  also  shows  the  minimum  bottom  angle  for 
rays  that  connect  with  the  sea  surface;  below  this,  data 
should  be  ignored.)  Interpreting  this  picture  now  as  a  bottom 
reflection  loss,  the  main  features  are  low  loss  up  to  a  critical 
angle  and  two  interference  lobes  beyond  this.  This  is  classic 
behavior  for  a  high  sound-speed  bottom  (see,  e.g.,  Jensen 
et  al. ,  1994  and  Nagl  et  al ,  1982);  once  the  critical  angle  is 
exceeded  there  is  the  possibility  of  simultaneous  reflections 
from  deeper  layer  boundaries  which  arrive  later  and  therefore 
interfere.  Indeed,  this  interference  pattern  is  evidence  that 
there  is  at  least  one  other  layer  boundary  present.  Further¬ 
more,  in  this  case  the  regularity  of  the  fringes  with  frequency 
suggests  exactly  two  boundaries.  The  drop  to  low  loss  at  the 
top  of  the  picture  is  an  artifact  caused  by  the  grating  lobe  that 
arcs  from  90°  at  1 500  Hz  (the  design  frequency)  to  40°  at  2 
kHz  in  Fig.  5(a).  Similarly,  the  low  losses  at  low  frequencies 
are  caused  by  the  inevitable  poor  angle  resolution  which  de¬ 
grades  the  up-to-down  ratio.  In  fact,  with  this  array  the 
beam  widths  to  3-dB  points  at  1500,  1000,  500,  200  Hz  are, 
respectively,  6.7,  10,  20,  50  degs.  in  retrospect,  all  these 
features,  particularly  the  two  fringes,  can  be  discerned  by  eye 
quite  clearly  in  the  original  array  response  plot. 

We  can  model  these  effects  by  applying  the  three-layer 
(two-boundary)  model  of  Sec.  II C  with  some  appropriate 
geoacoustic  parameters.  Geoacoustic  models  for  the  vicinity 
are  given  by  Max  et  al  (2000),  core  analysis  by  Tonarelli 
et  al.  (1993),  and  independent  inversion  results  from 
MAPEX2000bis  by  Siderius  et  al.  (to  be  published,  2002), 


As  shown  in  Table  II,  we  take  parameters  from  the  first  ref¬ 
erence  as  a  starting  point  and  then  make  minor  adjustments 
by  hand;  in  particular,  we  adjust  the  layer  thickness  for 
agreement  of  fringe  separation.  The  reflection  loss  is  shown 
in  Fig.  5(c).  If  we  then  simulate  the  effect  of  beamforming 
(see  the  Appendix),  the  result  is  remarkably  similar  to  the 
experimental  plot  [Fig.  5(d)]. 


B,  “Silt"  site:  N.  Elba,  1  December  2000 

Experimental  array  response  (AR)  and  reflection  loss 
(RL)  curves  for  a  silt  bottom  are  shown  in  Figs.  6(a)  and  (b). 
There  is  now  no  noise  notch  because  the  array  has  no  maxi¬ 
mum  sound  speed  above  it.  The  most  noticeable  feature  of 
the  reflection  plot  is  the  irregularity  of  the  interference 
fringes  and  the  features  at  low  grazing  angles.  In  fact,  we  see 
superimposed  fine  and  coarse  fringes  demonstrating  that 
there  must  be  more  than  two  layer  boundaries.  Furthermore, 
since  low  frequencies  tend  not  to  be  affected  by  thin  layers, 
we  can  see  that  the  lower  part  of  Fig.  6(b)  is  dominated  by 
fine  fringes  (200-Hz  separation)  caused  by  a  thick  layer  with 
critical  angle  —20°.  In  contrast,  the  coarser  fringes  seen 
starting  at  1200  Hz,  90°  must  be  caused  by  a  thinner  layer  (a 


FIG,  4.  Sound -speed  profiles  for  all  sites  derived  from  near  simultaneous 
CTDs  for  all  sites  except  New  Jersey,  which  was  XRT.  Labels  indicate  New 
Jersey  Shelf  (NiS),  Adventure  Rank  (Adv),  S,  Sicily  (Sic),  "sand"  site 
(solid),  “rock"  site  (dashed),  and  Elba  (Elh),  "silt"  site  (solid),  “silt"  site 
without  wind  (dashed),  “mud"  site  (dot-dash). 
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FIG,  5.  (a)  Experimental  VLA  beam  response  with  32  elements  #  0.5-m  separation  at  the  Ragusa  Ridge.  S.  Sicily,  ‘'sand”  site,  (b)  Reflection  loss  deduced 
from  the  up-io-down  ratio  of  the  beams  from  the  VLA  (min  bottom  angle  0°)*  (c)  Reflection  loss  from  a  simple  layer  model  with  parameters  adjusted  by  hand 
as  in  Table  II.  (d)  Heamibrmed  noise  resulting  from  the  modeled  reflection  loss  in  a  simple  sheet  noise  model  [Eqs.  (Al)t  (A2),  (A6)]. 


factor  of  3  or  4  thinner,  judging  by  the  fringes)  with  a  some* 
what  laiger  critical  angle.  This  already  provides  strong  clues 
to  the  bottom  structure.  Notice  lhat  the  evidence  for  multiple 
layers  is  only  manifest  above  about  900  Hz. 


In  earlier  experiments  at  nearby  sites,  Jensen  (1974) 
found  a  good  two-boundary  fit  at  low  frequencies  using 
normal-mode  propagation  modeling,  and  this  was  used  as  a 
starting  point  for  subsequent  inversion  work  by  Gingras  and 


TABLE  II.  Geoacoustic  parameters  for  all  sites  (initial  from  reference  and  final). 


Site 

Speed 

(m/s)1 

Thickness 

(my 

Density 

(gm/cc)a 

Attenuation 

(tiB/M* 

Ref* 

ct 

g2 

c3 

c4 

hi 

h3 

p2 

p3 

P4 

a2 

a? 

v4 

"silt" 

1690 

1750 

4.3 

,,, 

1.9 

... 

2.0 

0.4 

0.4 

TTMA 

Adv,  B. 

1511 

1630 

1650 

2.5 

... 

1*9 

-** 

1.9 

0.4 

0*4 

Fig  9 

“sand** 

... 

1554 

1950 

5 

... 

1,3 

2.0 

0.18 

0.3 

MFHTB 

Ragusa 

1512 

1554 

1800 

1.0 

,,, 

L2 

*•* 

2*0 

0.18 

♦  .. 

0.2 

Fig  5 

■'silt" 

... 

1600 

1530 

1600 

0.2 

2.5 

? 

1.75 

1,8 

7 

0.13 

0.15 

MWJ 

N.  Elba 

1514 

1600 

1530 

1600 

0,7 

2.1 

1.8 

1.75 

1.8 

0.13 

0*14 

0*15 

Fig  6 

"mud” 

1645 
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16(H) 

0.3 

6 

7 

1.5 

1.8 

7 

0*06 

0.15 

MWJ 

E.  Elba 

1509 

1530 

1471 

1530 

0,8 

3.5 

14 

1.2 

L8 

0*14 

0,06 

0.15 

Fig  7 

'■day" 

1530 

— 

? 

**• 

L4 

••• 

7 

0.1 

... 

7 

H 

NX  Sh. 

1499 

1695 

. . ■ 

1700 

3*0 

1.85 

•  *  * 

L95 

0,15 

* ., 

0,15 

Fig  10 

’Layers  are  numbered  1-4  from  the  top.  Layer  1  is  always  water,  but  only  the  sound  speed  is  shown.  Layer  4  is  always  a  half-space.  Sometimes  layer  3  is 
missing,  i.e.,  zero  thickness.  A  ‘T*  indicates  that  no  information  was  available  from  the  reference  cited. 

bBold  characters  indicate  values  actually  used  in  the  calculations  for  the  figure  indicated.  Normal  characters  indicate  values  suggested  by  the  references  cited. 
TTMA  =  Torture  Eli  it  aL  (1993)1  MFHTB  =  Max  etal.  (2 (XX));  MW J “Murphy  el  aL  (1976);  H=  Hamilton  (1980). 
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FIG.  6.  (a)  HxperimcntaJ  VLA  beam  response  with  32  elements  <3*  ().5-m 
separation  at  the  N.  Hlba  "silt”  site,  (b)  Reflection  loss  deduced  from  the 
up-to-down  ratio  of  the  beams  from  the  VLA  (min.  bottom  angle  5.09°).  (c) 
Bcamformcd  noise  resulting  from  the  reflection  loss  modeled  using  param¬ 
eters  adjusted  by  hand  as  in  Table  II.  (d)  Experimental  beam  response  at  the 
same  site  on  a  different  occasion  when  there  was  no  wind.  Note  extremely 
weak  "up”  and  "down”  contributions,  (c)  Anomalous  reflection  loss  de¬ 
duced  from  the  up-to-down  ratio  of  the  beams  from  the  VLA  (min.  bottom 
angle  6.03°).  The  ratio  tends  to  unity  since  the  up  and  down  steered  beams 
are  in  reality  dominated  by  sidelobes  of  the  strong  "horizontal”  noise  seen 
in  (d). 


Gerstoft  (1995).  Parameters  are  given  in  Table  II.  These  pa¬ 
rameters,  when  inserted  in  the  noise  model,  give  good  agree¬ 
ment  with  Fig.  6(b)  at  low  frequencies,  but  to  improve  agree¬ 
ment  at  high  frequencies  we  need  to  insert  a  thin  layer  with 
high  sound  speed  (Table  II).  In  fact,  there  is  a  lot  of  evidence 
from  core  analysis  for  such  a  thin,  high-velocity  layer  sitting 
on  top  of  the  sediment  near  this  site  [Akal  el  al.  (1972); 
Murphy  and  Olesen  (1974);  Jensen  (1974);  Murphy,  Wasilj- 


eff.  and  Jensen  (1976);  Holland  and  Osier  (2000)].  In  Fig. 
6(c)  we  show  a  very  good  fit  with  a  four-layer  model. 

C.  “Silt”  site,  no  wind:  N.  Elba,  29  November  2000 

The  silt  site  was  also  visited  on  an  occasion  when  the 
sea  was  flat,  calm,  and  there  was  hardly  any  wind  (^3  kt). 
There  is  a  dramatic  difference,  visible  immediately  from  the 
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FIG.  7,  (a)  Reflection  loss  deduced  from  the  up-to-down  ratio  of  the  beams 
from  the  VLA  al  the  E.  Elba  “mud"  site  (min,  bottom  angle  7.08").  (b) 
Beamformed  noise  resulting  from  the  reflection  kiss  modeled  using  param¬ 
eters  adjusted  by  hand  as  in  Table  IL 

AR  plot.  The  “up”  and  “down”  levels  are  understandably 
both  much  weaker,  but  there  are  still  sources  near  horizontal, 
probably  distant  shipping.  On  the  assumption  that  there  is  no 
vertical  noise  at  all,  the  beamformed  responses  are  really 
responses  of  sidelobes  to  this  horizontal  noise.  The  response 
is  therefore  virtually  the  same  up  as  down,  and  the  ratio 
tends  to  unity  (zero  dB),  This  is  essentially  what  we  see  in 
Figs.  6(d)  and  (e).  Although  there  is  structure  it  is  mainly 
anomalous.  However,  it  is  very  obvious  from  both  the  sym¬ 
metry  in  AR  plot  and  the  differences  in  level  between  verti¬ 
cal  and  horizontal  in  the  AR  plot  that  this  has  happened,  so 
the  technique  is  reasonably  robust  from  this  point  of  view.  It 
is  stressed  that  it  is  the  beamforming  that  is  breaking  down 
rather  than  the  theory  of  this  method:  the  up-to-down  ratio  is 
just  difficult  to  measure  in  the  presence  of  the  loud  horizon¬ 
tal  sources, 

D.  “Mud”  site:  E.  Elba,  30  November  2000 

The  mud  site  was  expected  to  be  lossy  with  sound 
speeds  lower  than  that  in  water  [Murphy,  WasiljefT,  and 


FIG,  8,  Reflection  loss  deduced  from  the  up-to-down  ratio  of  the  beams 
from  the  VLA  at  the  Ragusa  Ridge,  S.  Sicily,  “rock"  site  (min.  bottom  angle 
2.47°). 

Jensen  (1976)],  In  Fig,  7{a)  there  is  evidence  of  higher  losses 
at  low  angles,  i.e.,  an  absence  of  a  clear  critical  angle.  Again, 
there  are  two  superimposed  fringe  patterns,  one  fine — with 
seven  or  eight  loss  maxima  visible  at,  say  90°,  and  one 
coarse — with  only  one  maximum  and  one  minimum  visible. 
This  suggests  three  boundaries,  two  of  which  are  wide  apart 
and  two  of  which  are  dose  together,  Murphy,  Wasiljeff,  and 
Jensen  (1976)  demonstrate  that  there  is  again  a  thin,  high- 
velocity  layer  at  this  site  but  now  in  the  middle  of  the  sedi¬ 
ment  layer.  The  result  of  a  hand  search  with  a  three- boundary 
model,  assuming  a  thin  layer  above  the  slow  sediment  layer, 
including  beamforming,  is  shown  in  Fig.  7(b)  (parameters 
are  in  Table  II),  By  setting  first  sound  speeds,  then  layer 
thicknesses,  then  density  and  absorption,  it  is  possible  inde¬ 
pendently  to  control  the  two  sets  of  fringes,  their  relative 
amplitudes,  and  the  low  angle  dependence.  The  search  can 
be  extended  to  a  fourth  boundary,  simulating  the  position  of 
the  thin,  fast  layer  within  the  main  sediment  layer,  but  this 
does  not  appear  to  improve  the  agreement. 

E.  “Rock"  site:  S,  Sicily  (Ragusa  Ridge),  23 
November  2000 

This  site  is  included  here  for  completeness,  but  we  sus¬ 
pect  that,  as  in  Figs.  6(d)  and  (e),  the  wind  was  not  strong 
enough  to  avoid  arti facts,  being  5-6  kts.  Nevertheless,  in 
Fig.  8  one  can  still  see  weak  structure  resembling  that  seen  in 
earlier  pictures  and  evidence  of  interference  through  layer¬ 
ing,  for  instance  an  interrupted  thin  layer  of  silt  above  rock. 
Although  one  would  expect  rock  to  be  a  strong  reflector,  a 
half-space  with  impedance  twice  that  of  water  would  still 
have  a  reflection  loss  of  9  dB  or  so  at  high  angles. 

F.  “Silt”  site:  SW  Sicily  (Adventure  Bank),  6  May  1999 

The  peak  values  of  reflection  loss  at  this  site  [Fig.  9(a)] 
appear  to  be  weaker  than  at  the  earlier  ones.  Suitable  envi¬ 
ronmental  data  are  provided  by  Tonarelli  et  al  (1993);  Caiti 
et  ai  (1996);  Caiti  (1996);  Snellen  et  ai  (2001);  Siderius 
et  ai  (2001)  (see  Table  II),  There  is  quite  a  spread  of  param- 
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FIG.  9.  (a)  Reflection  loss  deduced  from  the  up-to-down  ratio  of  the  beams 
from  the  VLA  at  the  Adventure  Bank.  S.  W  Sicily,  “sill”  site  (zero  bottom 
angle  corresponds  to  2.80”  at  the  array),  (b)  Beamformed  noise  resulting 
from  the  reflection  loss  modeled  using  parameters  adjusted  by  hand  as  in 
Table  II. 


FIG.  10.  (a)  Reflection  loss  deduced  from  the  up-lo-down  ratio  of  the  beams 
from  the  VLA  at  the  New  Jersey  Shelf  “stiff  clay”  site  (min  bottom  angle 
0°).  (b)  Beamformed  noise  resulting  from  the  reflection  loss  modeled  using 
parameters  adjusted  by  hand  as  in  Table  II. 


eters  between  these  references,  partly  because  of  real  spatial 
variability,  with  sediment  densities  from  1.6  and  1.9  and 
sound  speeds  from  1670  and  1750  in  the  first  few  meters. 
The  model  fits  shown  in  Fig.  9(b)  are  good  but  not  quite  as 
good  as  in  the  other  cases.  It  is  possible  that  deeper  high¬ 
speed  layers  could  account  for  the  “lamb  chop”  shape  of  the 
experimental  fringes  (by  providing  loss  peaks  at  the  higher 
critical  angle),  but  in  retrospect  these  shapes  look  rather  sus¬ 
piciously  unlike  theoretical  curves.  The  prime  suspect  is  the 
addition  of  nonacoustic  noise  which  artificially  reduces  the 
up-to-down  ratio.  It  is  known  that  these  data  were  affected 
by  a  form  of  uncorrelated  (pink)  noise,  namely,  quantization 
noise  when  the  gains  had  been  set  for  experimental  purposes 
other  than  gathering  noise. 

An  interesting  separate  point,  taken  up  in  Harrison  and 
Simons  (2001),  is  that,  although  there  are  changes  over  sev¬ 
eral  hours  in  noise  directionality  that  are  obvious  in  the  beam 
response,  the  derived  reflection  loss  is  more  or  less  immune 
to  them. 


G.  “Stiff  clay”  site:  New  Jersey  Shelf,  17  May  2001 

There  appears  to  be  little  quantitative  geoacoustic  data 
available  for  this  site,  although  it  is  thought  to  be  a  stiff  clay 
(Goff  et  al.y  1999).  Hamilton  (1980,  1987)  provides  some 
average  properties  for  “continental  terrace  (shelf  and  slope) 
environment.”  Between  clayey  silt  and  silty  clay,  sound 
speeds  are  1546-1517  m/s  and  densities  are  1.489-1.480 
g/cc.  Silt  clays  have  a  low  loss  of  about  0.1  dB/wavelength. 
Initial  parameters  are  shown  in  Table  II. 

The  experimental  reflection  loss  curve  in  Fig.  10(a) 
shows  some  layer  structure,  but  the  main  feature  is  a  simple 
step  at  around  35°  up  to  about  10-dB  loss.  This  can  be  mod¬ 
eled  either  as  a  half-space  or  as  a  multiple  boundary.  The 
sound  speed  in  the  basement  is  determined  by  the  obvious 
critical  angle,  and  the  layer  thickness  is  determined  by  the 
fringe  separation.  To  get  a  good  fit  to  the  depth  of  modula¬ 
tion,  while  treating  the  layer  as  clay,  it  was  necessary  to 
choose  a  density  close  to  that  of  water  (1.05),  but  this  con¬ 
flicts  with  the  assumption  of  a  clay  sediment.  Equally  con¬ 
vincing  fits  can  be  obtained  by  assuming  either  a  half-space 
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of  fine  sand  with  sound  speed  1700  m/s  and  density  1,8  g/ee 
(Hamilton,  1987),  or  a  3-m  layer  and  basement  of  slightly 
differing  properties  close  to  that  of  fine  sand.  The  effect  of 
beamforming  with  the  parameters  of  Table  II  is  shown  in  Fig. 
10{b).  Note  that  this  example  happens  to  reveal  an  ambiguity 
or  nonuniqueness  in  the  inversion  process  which  is  common 
to  many  inversion  methods.  The  experimental  noise  tech¬ 
nique  deduces  a  reflection  loss  which  may  indeed  have  mea¬ 
surement  errors,  but  the  reflection  properties  are  not  ambigu¬ 
ous  in  the  same  sense  as  the  geoacoustic  properties,  A 
sensitivity  analysis,  particularly  in  this  case,  would  be  inter¬ 
esting  but  it  is  felt  to  be  outside  the  scope  of  this  paper.  Some 
of  these  issues  are  discussed  in  Sec.  V. 

H.  Data  quality 

The  dynamic  range  of  the  reflection  loss  {i.e.,  depth  of 
modulation  of  fringes)  sometimes  appears  to  be  slightly 
smaller  in  the  measurement  than  in  the  calculation.  An  obvi¬ 
ous  possible  cause  is  the  addition  of  {frequency -depen dent) 
uncone lated  noise  which  potentially  spoils  the  up-to-down 
ratio.  It  has  the  effect  of  disproportionately  boosting  the 
weak,  upward -going  signal,  thus  making  the  bottom  appear 
to  be  a  stronger  reflector  In  terms  of  the  cross-spectral  den¬ 
sity  matrix  for  the  array,  uncorrelated  noise  adds  to  the  diag¬ 
onal  only.  Adding  the  equivalent  of  10^  of  the  downward 
power  to  both  up  and  down  would  result  in  a  10-dB  ceiling 
in  measured  RL. 

A  form  of  unconelated  noise  (in  the  ADVENT99  data 
case  only)  was  quantization  noise.  Another  contender  is  rat¬ 
tling  or  knocking  of  hydrophones  in  their  oil-filled  pipe.  The 
latter  was  found  on  several  occasions  and  it  has  been  pos¬ 
sible  to  find*  at  one  site  (New  Jersey  Shelf)*  a  morning  with 
intermittent  knocking  [multichannel  time  series  in  Fig,  11(a)] 
and  an  afternoon  without  [Fig.  1 1  (b)].  The  surprising  finding 
is  that,  although  the  noise  directionality  has  changed  between 
morning  and  afternoon,  this  degree  of  interference  makes  no 
visible  difference  to  the  corresponding  derived  reflection 
loss;  Fig.  12  shows  the  equivalent  of  Fig.  10(a)  but  with 
knocking.  Evidently,  although  the  raw  noise  time  series  looks 
severely  contaminated  in  the  morning*  the  added  absolute 
power  is  insignificant. 

V.  DISCUSSION  OF  PERFORMANCE 
A*  Strengths 

The  experimental  data  demonstrate  a  number  of 
strengths  of  the  noise  inversion  method. 

(i)  No  sound  source  is  required, 

(ii)  No  propagation  or  noise  model  is  required, 

(iii)  The  method  tolerates  unknown  source  distributions 
and  noise  directionality.  It  is  self-compensating. 

(iv)  The  method  tolerates  distant  point  sources  and  ships. 

(v)  The  method  tolerates  refracting  environments. 

(vi)  The  method  tolerates  range-dependent  environments, 

(vii)  In  contrast  with  most  other  inversion  methods,  high- 
frequency  performance  is  good, 

(viii)  The  method  explicitly  shows  the  critical  angle  and  an 
indication  of  the  number  of  layers;  one  can  distin- 


FIG.  ii.  (a)  Complete  64  channels  of  lime  series  (AM,  1 7-05 -01)  showing 
jxHtr  quality  data  with  sound  emanating  from  points  on  the  array.  Vertical 
channel  separation  is  proportional  to  hydrophone  separation,  (b)  Complete 
64  channels  of  lime  series  {PM,  174)5-01)  showing  good  quality  daia 


guish  by  eye  between  one,  two,  and  more  than  two 
bottom  layers.  Therefore,  it  is  possible  to  isolate  by 
hand  the  sound  speed  and  thickness  in  the  top  few 
layers.  The  depth  of  modulation  is  controlled  by  the 
relative  strengths  of  the  boundary  reflections,  i.e., 
their  impedance  mismatches,  so,  having  fixed  the 
sound  speeds  one  can  determine  the  densities. 

(ix)  The  method  is  tolerant  to  poor  data  quality. 

An  interesting  question  is,  if  the  bottom  loss  varies  with 
position,  does  this  method  measure  a  spatial  average  bottom 
loss  or  is  it  spatially  biased?  Ray  theory  (see  the  Appendix) 
suggests  that  the  crucial  bottom  area  is  quite  near  to  the 
receiver.  In  fact,  it  should  be  in  the  vicinity  of  the  first  bot¬ 
tom  bounce  which,  of  course,  is  a  function  of  angle.  This 
hypothesis  has  been  confirmed  by  conducting  a  numerical 
experiment  with  a  range-dependent  wave  model  in  which  a 
variable  size  area  beneath  the  array  has  distinguishable  re¬ 
flection  properties  from  more  distant  points  (Harrison  and 
BaJdacci,  2002).  The  corollary  is  that  a  drifting  array  could 
use  this  method  to  survey  spatial  bottom  variations.  Indeed, 
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FIG.  12.  Reflection  loss  deduced  from  the  up-to-down  ratio  of  the  beams 
from  the  VLA  at  the  New  Jersey  Shelf  “stiff  clay”  site  but  using  poor 
quality  data  as  shown  in  Fig.  1 1  (a).  Note  that  dilTerences  between  this  and 
the  good  quality  data  result  [Fig.  10(a)]  are  barely  perceptible. 


ior  of  the  physical  parameters.  An  additional  weakness  in  the 
noise  case  is  that  we  have  lost  all  absolute  phase  informa¬ 
tion,  although  we  have  retained  some  relative  phase  informa¬ 
tion  in  the  form  of  the  interference  fringes.  If  we  were  to 
Fourier  transform  a  complex  (frequency-dependent)  reflec¬ 
tion  coefficient  we  would,  in  principle,  obtain  an  impulse 
response  of  the  layers  akin  to  standard  geophysical  displays 
[see  the  “move-out"  technique  in  Holland  and  Osier  (2000)]. 
The  inversion  ambiguity  would  therefore  be  the  same  as  for 
those  techniques.  If,  on  the  other  hand,  we  Fourier  transform 
the  modulus  square  of  the  reflection  coefficient  we  would 
obtain  the  autocorrelation  function  of  the  impulse  response. 
Thus,  layer  spacing  information  is  still  there  but  it  is  mixed 
up  to  the  extent  that  one  cannot  distinguish  the  order  of,  say, 
a  thick  layer  and  a  thin  layer.  Nevertheless,  even  with  a 
power  reflection  coefficient,  the  evident  critical  angle,  inten¬ 
sity  contrast,  and  fringe  separation  mean  that  sound  speed, 
impedance  mismatch,  and  layer  thicknesses  are  fairly  well 
determined  as  long  as  we  assume  a  reasonably  small  number 
of  layers. 


this  drift  technique  has  been  tried  in  a  recent  experiment 
(BOUNDARY2002,  April  2002,  unpublished),  again  using 
the  16-m  central  part  of  the  VLA.  Convincing  bottom 
changes  were  seen  during  a  several  mile  drift. 

B.  Analysis  issues 

1.  Random  errors 

Assessing  the  impact  of  errors  quantitatively  is  a  diffi¬ 
cult  problem,  and  the  purpose  of  Harrison  and  Baldacci 
(2002)  was  to  investigate  some  of  these  problems  by  numeri¬ 
cal  simulation.  Among  the  problems  discussed  elsewhere  in 
this  paper  are:  addition  of  uncorrelated  noise  to  the  hydro¬ 
phones;  deviation  of  the  array  from  the  vertical;  sidelobe 
leakage  in  the  presence  of  powerful  directional  sources;  er¬ 
rors  in  up-down  ratio  for  small  angles  where  low  losses  are 
expected;  resolvable  eigenrays  from  nearby  shipping. 

2.  Systematic  errors 

Two  effects  are  obvious  in  all  the  comparisons  between 
experimental  and  modeled  reflection  loss  (including  the 
OASR,  OASN  pair  in  Fig.  2).  One  is  that,  inevitably,  at  angles 
as  small  as  the  beam  resolution  the  up-to-down  ratio  tends  to 
unity.  Thus,  there  is  a  predictable  region  near  horizontal 
where  the  result  is  anomalous.  The  region  gets  wider  as  fre¬ 
quency  lowers,  but  improvements  are  possible  either  by  us¬ 
ing  techniques  such  as  adaptive  beamforming  or  by  increas¬ 
ing  the  aperture.  Another  effect  that  depends  on  hydrophone 
separation  rather  than  array  size  is  the  onset  of  a  grating  lobe 
at  the  design  frequency.  This  typically  results  in  anomalously 
low  reflection  loss. 


C.  Experimental  issues 

In  general,  there  are  a  number  of  experimental  issues. 

1.  Requirement  for  a  VLA 

It  is  true  that  the  existing  experimental  arrangement  us¬ 
ing  a  62-meter  VLA  is  cumbersome  and  expensive.  How¬ 
ever,  the  aperture  used  in  all  the  cases  shown  here  is  only 
15.5  meters  (32  elements  @  0.5-m  separation),  and  one 
could  contemplate  two  alternatives.  One  is  an  expendable 
VLA  that  is  left  to  drift  in  much  the  same  way  as  a  sonobuoy. 
The  other  is  a  drifting  or  moored  synthetic  aperture,  where 
one  hydrophone  is  fixed  in  depth  and  the  other  moves  up  and 
down.  In  this  way  one  could  obtain  vertical  coherence,  or 
rather,  cross-spectral  density  for  all  necessary  separations  as 
required  by  Eq.  (9). 

2.  Up- down  beamforming  contrast 

Although  the  interference  patterns  seen  in  the  experi¬ 
mental  curves  resemble  those  seen  in  theory,  the  theoretical 
reflection  loss  curves  often  contain  taller  spikes.  This  effect 
is,  of  course,  highlighted  by  displaying  in  dBs,  but  there  are 
various  possible  causes  that  could,  in  principle,  be  serious. 
Possibilities  include  addition  of  nonacoustic  noise  (electrical, 
quantization)  or  uncorrelated  noise  (strumming,  knocking)  or 
poor  beamforming  (incorrect  assumptions  on  hydrophone 
gain,  hydrophone  separation  or  location,  poor  angle  discrimi¬ 
nation,  poor  knowledge  of  sound-speed  profile).  Most  can  be 
checked  by  inspection  of  co-K  plots,  intensity  histograms  of 
individual  channels,  or  the  cross-spectral  density  matrix  near 
the  diagonal. 


3.  Inversion  problems 

Even  given  an  error-free  measurement  of  (power)  reflec¬ 
tion  loss,  the  geoacoustic  properties  that  result  from  inver¬ 
sion  are  not  necessarily  unique.  This  is  common  and  well 
known  with  many  inversion  techniques,  and  the  usual  solu¬ 
tion  is  to  introduce  constraints  on  the  magnitudes  and  behav- 


3.  Array  tilt 

Even  if  the  tilt  is  known,  it  is  difficult  to  compensate 
mathematically  without  knowing  the  horizontal  noise  source 
distribution — but  to  know  this  defeats  the  object  of  this 
method,  i.e.,  simplicity.  Therefore,  we  need  to  make  sure  that 
the  tilt  stays  inside  some  bounds  set  by  the  array’s  angle 
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resolution.  However,  fortuitously,  when  the  resolution  is 
poor  {for  low  frequencies  and  steep  angles)  the  reflection 
loss  varies  slowly  with  angle.  This  is  no  coincidence;  the 
former  depends  on  the  vertical  scale  of  the  may,  whereas  the 
latter  depends  on  the  vertical  scale  of  the  layers.  In  fact,  a 
glance  at  any  of  the  reflection  loss  plots  shows  relatively 
gradual  changes  with  angle,  so  that  a  tilt-induced  spread  of  a 
few  degrees  cannot  produce  significant  effects  (Harrison  and 
Baldacci,  2002), 


thanks  are  due  to  Martin  Siderius  and  Charles  Holland  for 
help  in  choosing  sites.  The  first  author  acknowledges  a  con¬ 
versation  with  David  Chapman,  Defense  Research  Establish¬ 
ment  Atlantic  (DREA),  on  this  technique  long  before  having 
the  opportunity  to  do  experiments, 

APPENDIX:  REFLECTION  LOSS  RELATION  TO  NOISE 
DIRECTIONALITY  IN  A  RANGE-DEPENDENT 
ENVIRONMENT 


4,  Low-angle  (near-horizontal)  data 

In  this  approach  the  reflection  loss  for  very  small  angles 
depends  on  small  differences  between  large  quantities.  It  is 
therefore  intrinsically  unreliable.  The  true  reflection  loss 
must  tend  to  zero  at  0°  so  one  could,  perhaps,  make  use  of 
this  fact  and  loss  values  just  inside  the  critical  angle  to  de¬ 
duce  the  absorption  of  the  top  sediment  layer, 

VI  CONCLUSIONS 

A  simple  technique  for  extracting  geoacoustic  param¬ 
eters  from  ambient  noise  measurements  has  been  proposed 
and  demonstrated  with  experimental  data.  Unlike  coherence- 
based  methods  that  use  a  hydrophone  pair,  this  technique 
does  not  require  a  noise  model.  It  can  therefore  still  work 
when  the  noise  directionality  and  bathymetry  are  unknown. 
The  workable  frequency  range  is  determined  by  the  array 
design  frequency  at  high  frequency,  where  grating  lobes 
spoil  the  up-to-down  beam  ratio,  and  by  its  angle  resolution 
at  low  frequency.  In  practical  systems  this  complements  con¬ 
ventional  active  inversion  techniques  where  sound- speed 
fluctuations  impede  model  matching  at  high  frequencies 
(Siderius  et  ct/.,  2001), 

A  theoretical  justification  was  given  in  terms  of  acoustic 
flux  and  also  ray  theory.  This  was  supported  by  the  results  of 
a  separate  numerical  study,  using  the  wave  models  OASN  and 
oasr,  in  which  there  was  close  agreement  between  the  up- 
to-down  beam  ratio  and  the  hypothesised  “true”  plane  wave 
bottom  loss, 

A  major  point  of  the  paper  is  that  there  are  very  clear 
variations  from  site  to  site  that  are  unaffected  by  temporal 
changes  in  noise  directionalities  at  those  sites.  In  addition, 
the  reflection  loss  interference  patterns  can  be  simply  mod¬ 
eled  with  parameters  close  to  those  already  found  in  other 
studies  for  the  same  sites.  So  far,  the  hand-searching  tech¬ 
nique  has  relied  on  inspecting  the  experimental  plot  for  criti¬ 
cal  angle,  fringe  spacing,  and  depth  of  modulation  to  find  a 
solution.  It  is  felt  that  the  method  is  now  ripe  for  an  auto¬ 
mated  search. 

The  method  appears  to  be  robust  against  most  common 
forms  of  experimental  data  contamination  except,  of  course, 
for  absence  of  wind, 
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Here,  we  choose  a  ray  approach  in  order  to  allow  for 
spatially  varying  environments  and  source  distributions  as 
already  been  developed  by  Harrison  (1997a).  Given  a  noise 
directionality  N(  <£,  S)  (a  function  of  elevation  angle  0  and 
azimuth  <fi)  and  the  vertical  array’s  beam  pattern  B(6,0Q)  for 
each  steer  angle  0o  t  the  array  response  Ai  0o)  is  given  by 


A(0O)  =  J  J  (A I) 

In  a  range-dependent  environment  the  directionality  jV(<£,  0) 
can  be  expressed  as  two  multiplicative  terms  [Harrison 
(1997a)];  in  the  equivalent  formula  for  coherence  the  beam 
pattern  is  simply  replaced  by  a  phase  term 

N(M  =  Q(e)F(<t>,0).  (A2) 

The  function  F  is  generally  a  complicated  function  of 
bathymetry,  sound  speed,  and  noise  source  distribution  in¬ 
cluding  multiple  arrivals  for  each  angle,  although  under  cer¬ 
tain  conditions  it  can  still  be  expressed  in  closed  form  [Har¬ 
rison  (1997a)].  For  instance,  it  is  a  geometric  series  in  a 
range -independent  environment.  The  function  Q  is  simply 

Q(8)=exp(-asp),  0^0 

~Rh(0h)exp(  —a(sc—sp))f  0< 0,  (A3) 

where  Ry{B\A  the  local  (within  one  ray  cycle  of  the  re¬ 
ceiver)  bottom  (power)  reflection  coefficient,  and  Of,  is  re¬ 
lated  by  Snell's  law  to  0  at  the  receiver.  The  exponential 
terms  represent,  respectively,  the  attenuations  along  the  re¬ 
sidual  parts  of  the  ray  directly  from  the  surface  (sp)  and 
directly  from  the  bottom  (sc  —  sp).  If  we  could  access  the 
directionality  itself,  dividing  the  downward  N  by  the  upward 
N,  we  would  eliminate  the  function  F,  leaving 


N(4>,-B)  Q(~8) 

mJTff)  =  Qi TV) =  Rb(  )e*p(  - ' -  2s»)l 

(A4) 


Since  this  is  true  for  all  azimuths  0,  it  is  also  true  for  the 
azimuth  integral  Assuming  that  the  effects  of 

absorption  over  one  ray  cycle  are  small,  we  obtain 


fAT(0.-W_ 


(A5) 


Given  appropriate  angular  resolution  in  Eq.  (Al),  we  can 
actually  measure  as  a  function  of  ft  and  hence  0h ,  since 
the  beam  response  tends  to  the  noise  directionality  N. 

From  the  ray  point  of  view  the  reason  for  this  simplicity, 
despite  inclusion  of  spatial  variation  of  the  environment  and 
noise  source  distribution,  can  be  seen  as  follows.  Imagine  a 
distant  point  noise  source  from  which  a  single  ray  extends. 
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after  multiple  surface  and  bottom  reflections,  to  the  receiver. 
For  every  such  ray  arriving  directly  (i.e.,  most  recently)  from 
the  surface,  there  is  a  corresponding  ray  arriving  from  the 
seabed  with  exactly  one  extra  bottom  reflection,  the  same 
number  of  surface  reflections,  and  almost  the  same  horizon¬ 
tal  wave  number.  Provided  the  environment  changes  only 
slowly  over  one  ray  cycle,  the  slight  shifts  in  position  of 
reflection  points  have  no  effect.  So,  all  distant  point  sources, 
regardless  of  their  strength,  will  appear  weaker  by  Rf>{^b) 
for  downward  tfthan  for  upward  0.  In  addition,  we  can  relax 
the  '"distant  source"  condition  provided  the  source  distribu¬ 
tion  is  close  to  uniform  as  in  the  subsequent  part  of  Harrison 
(1997a).  The  reasoning  here  is  that  with  a  sheet  source  there 
is  no  geometric  spreading  (it  cancels  out),  so  the  only  differ¬ 
ence  between  an  upgoing  and  a  downgoing  ray  at  exactly  the 
same  angle  is  the  nearby  bottom  reflection.  The  same  argu¬ 
ment  applies  even  if  there  have  been  several  ray  cycles  be¬ 
tween  the  surface  source  and  the  receiver. 

In  order  to  simulate  the  effects  of  beamforming  on  the 
directionality  [see  Eq.  (Al)]  we  note  that,  in  the  case  where 
the  noise  sources  are  uniformly  distributed  dipoles  and  the 
environment  is  independent  of  range,  the  function  F  in  Eq, 
(A2)  reduces  to 

F(  0, 0)  =  sin  0,  f[  1  -  Rt(  0x)Rh(  0*)exp(  -  as,)  ]t  (A6) 

with  R s  and  being  surface  loss  and  surface  angle,  respec¬ 
tively.  We  can  therefore  estimate  the  array  response  using 
Eqs.  (Al)  and  (A2)  knowing  the  beam  pattern  (from  the 
array  dimensions).  Dividing  downward  by  upward  array  re¬ 
sponse,  we  obtain  the  ratio  of  beam  responses  which  is  an 
approximation  to 

A{-0O) 

<A7> 

We  note  that,  although  F  does  not  truly  cancel  out,  its  influ¬ 
ence  on  the  ratio  9\(0h(0o))  is  weak.  Therefore,  we  assume 
that  an  estimate  of  F  under  the  ideal  conditions  of  a  uniform 
noise  source  distribution,  as  in  Eq.  (A6),  will  suffice, 
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To  invert  for  parameters  of  the  seabed  and  the  w  ater  column,  global  optimisation  methods 
are  used.  These  methods  search  for  the  best  parameter  combination,  so  as  to  make  the 
difference  between  modelled  and  measured  pressure  fields  as  small  as  possible.  Simulated 
annealing  and  genetic  algorithms  (GA)  are  commonly  used  as  global  optimisation  algorithms 
for  this  inversion.  Their  efficiency  depends  on  the  tuning  of  the  algorithm.  For  a  bottom 
model  with  in  the  order  of  ten  geo-acoustic  parameters,  a  relatively  new  global  optimisation 
method  -  Differential  Evolution  (DE)  -  is  applied.  We  claim  a  ten  times  higher  efficiency  of 
DE.  as  it  finds  the  optimal  parameter  combination  with  only  one  tenth  of  the  number  of  model 
evaluations  needed  by  GA.  Apart  from  the  tuning  of  the  algorithms,  the  performance  of  DE  in 
comparison  with  GA  is  investigated  in  terms  of  efficiency  and  robustness. 


1.  INTRODUCTION 

In  the  last  two  decades  global  optimisation  methods,  such  as  simulated  annealing  (SA)  and 
genetic  algorithms  (GA),  have  been  applied  to  solve  inversion  problems  in  underwater 
acoustics.  In  this  inversion  a  pressure  field  calculated  by  a  model  is  compared  with  a 
measured  pressure  field,  using  a  cost  or  energy  function.  The  search  is  for  that  parameter 
combination  that  maximises  similarity  between  modelled  and  measured  field,  i.e.,  giving  the 
lowest  energy  value.  In  general,  a  fairly  large  amount  of  parameters  arc  sought,  typically  5  to 
20.  This  is  called  ‘matched  field  inversion’. 

A  purely  random  search  would  require  too  many  forward  model  calculations.  Search 
methods  strongly  reduce  this  number  of  function  evaluations,  but  have  the  risk  to  fail  by 
being  trapped  in  a  local  minimum  of  the  energy  function.  To  overcome  this  problem,  global 
optimisation  methods,  such  as  simulated  annealing  and  genetic  algorithms,  have  been 
designed. 


Performance  criteria  for  optimisation  methods  are  efficiency,  robustness  and  accurateness. 
Since  the  forward  model  calculations  take  most  of  the  computation  time  during  optimisation, 
efficiency  is  measured  by  the  required  amount  of  these  calculations,  typically  a  few  tens  of 
thousands.  The  stochastic  nature  of  these  methods  causes  that  even  a  global  search  run  can 
get  stuck  in  a  local  minimum.  Robustness  is  defined  as  the  probability  of  an  optimisation  run 
to  end  in  the  global  minimum,  i.e.,  the  probability  of  success  of  the  optimisation  method. 
Accurateness  has  to  do  with  the  distance  between  the  outcome  of  a  successful  optimisation 
run  and  the  precise  location  of  the  global  minimum.  In  this  article  the  attention  is  focused  on 
robustness  and  efficiency. 

Optimisation  methods  have  a  few  setting  parameters  (e.g.  crossover  probability  and 
mutation  rate  in  GA),  that  have  to  be  optimised  themselves.  It  will  be  shown  that  these 
settings  strongly  influence  the  performance  of  the  method.  Although  the  optimal  setting  can 
be  problem  specific,  we  hope  to  find  optimum  settings  that  only  need  minor  adjustments  for 
other  problems. 

It  is  already  well  established  that  in  matched  field  inversion  GA  outperforms  SA. 
However,  there  is  a  need  for  methods  even  better  than  GA.  A  relatively  new  global  search 
method,  differential  evolution  (DE),  is  introduced.  Its  performance,  in  terms  of  efficiency  and 
robustness,  is  shown  to  be  better  than  that  of  genetic  algorithms. 


2.  BASIC  PRINCIPLES  OF  GENETIC  ALGORITHMS  AND  DIFFERENTIAL 

EVOLUTION 

In  underwater  acoustics  the  GA  has  been  introduced  by  Gerstoft  [1].  Its  principle  of 
operation  is  summarised  as  follows.  To  optimise  an  energy  function  E,  GA  improves 
populations  of  q  parameter  combinations  during  successive  generations.  A  GA  optimisation 
run  starts  by  creating  at  random  an  initial  population  of  q  members.  Each  member  represents 
a  certain  parameter  value  combination  m  and  its  energy  E  is  calculated.  A  population 
member  with  a  lower  E  (a  higher  fitness)  has  a  higher  probability  to  be  selected  in  one  (or 
more)  of  the  pairs  of  parents  to  be  formed.  This  results  in  q! 2  pairs  of  parents  with  a  larger 
proportion  of  fit  members.  From  each  pair  of  parents  two  children  arc  created  by  a  random 
exchange  of  parental  values  (steered  by  the  crossover  probability  pc),  followed  by  random 
changes  of  individual  parameter  values  (with  mutation  probability  pm).  A  set  of  q  children 
vectors  results.  The  next  generation  is  established  by  taking  at  random  /rg  members  of  these 
children  and  selecting  the  (1  -  fr)q  most  fit  members  of  the  current  population.  This  process 

is  repeated  over  several  hundreds  of  generations. 

GA’s  have  the  following  setting  parameters:  population  size  q,  crossover  probability  pc, 
mutation  probability  p„„  reproduction  size  fr  and  the  number  of  generations  No- 

Differential  evolution,  just  like  GA,  starts  with  an  initial  population  of  q  randomly  chosen 
parameter  value  combinations  m.  These  m's  arc  improved  during  successive  generations  of 
constant  size  q,  in  the  sense  that  a  descendant  replaces  an  m,  becoming  its  successor,  if  it  has 
a  lower  energy.  The  distinctive  feature  of  DE  is  the  way  in  which  these  descendants  are 
created.  In  [2]  various  ways  to  generate  new  m’s  are  described.  Here  only  the  following 
procedure  is  considered.  At  the  start  of  generation  k  the  parameter  vectors  mA ,,...,  mA?  arc 

given  and  for  each  of  them  a  descendant  is  created,  being  a  potential  successor.  To  create  this 
descendant  dA .  of  mt  j  first  a  partner  pt  (  is  constructed  according  to 


?k,,=mkh+F(mklj~mkh) 

with  three  different  vectors  and  chosen  at  random  from  the  population  and  F 

a  scalar  multiplication  factor  between  0  and  1.  The  descendant  <  of  m,,  results  from 
applying  crossover  to  mA .  and  p*  (  with  crossover  probability  pc.  A  higher  value  of  pc  leads 
(on  the  average)  to  more  dimensions  of  pkJ  to  be  copied  into  dt(,  while  the  values  for  the 
remaining  dimensions  arc  taken  from  mir  dt,  only  replaces  mt/,  becoming  its  successor, 
if  its  energy  is  lower. 

The  setting  parameters  of  DE  are  population  size  q,  multiplication  factor  F,  crossover 
probability  pc  and  the  number  of  generations  Nc- 

Comparing  GA  and  DE,  both  show  improvement  of  a  fixed  sized  population  of  parameter 
combinations.  However,  in  each  generation  of  GA  a  fixed  fraction  of  the  population  is 
replaced  by  children,  whether  they  arc  better  or  not.  In  DE  each  member  of  the  population 
will  be  replaced  by  its  descendant  only  if  this  descendant  has  lower  energy.  On  the  other 
hand,  the  fitness  of  a  member  increases  the  probability  of  that  member  to  play  a  role  in  the 
creation  of  children  in  GA,  where  in  DE  this  fitness  docs  not  influence  the  probability  of  a 
member  to  be  selected.  Another  difference  is  the  way  a  partner  pijt  is  found  for  m* , .  While 

in  GA  the  partner  is  an  element  of  the  given  generation,  in  DE  it  is  constructed  from  the 
population  and  thus  not  an  element  of  this  generation.  The  partner  creation  process  in  DE  is 
steered  by  the  multiplication  factor  F  and  replaces  the  mutation  mechanism  of  GA. 


3.  PERFORMANCE  OF  GENETIC  ALGORITHMS  AND  DIFFERENTIAL 

EVOLUTION 

A  test  function  given  in  [3]  is  known  to  have  characteristics  that  correspond  to  those  of  the 
geo-acoustic  inversion  problem.  Its  features  are:  number  of  parameters  Np  =  6,  various  local 
minima,  parameter  coupling  and  a  varying  sensitivity  to  the  parameters  /«<.  The  global 
minimum  is  at  the  origin  m,  =  0.  The  research  on  the  performance  of  GA  and  DE  presented 
here,  is  restricted  to  minimisation  of  this  energy  function  called  £/. 

To  determine  the  probability  of  success  multiple  optimisation  runs  arc  needed,  each  run 
starting  with  a  different  randomly  chosen  population  of  q  m's.  The  number  of  independent 
runs  per  setting  is  denoted  Nr.  Nc  is  the  amount  of  these  runs  that  have  converged  to  the 
global  minimum.  The  fraction  of  converged  runs  serves  as  an  estimate  for  the  probability  of 
success  p„  according  to 

Ps=lt  ^  (Tp‘=~N^  an<1  =-JP^~Ps)Nr 

When  comparing  the  performance  of  different  settings  of  the  optimisation  algorithms,  also 
the  statistical  uncertainty  on  ps,  being  crPi »  should  be  taken  into  account  (where  N(  has  a 

binomial  distribution).  The  efficiency  is  measured  by  the  number  of  function  evaluations  or 
forward  model  calculations  N,  given  by  N  =frq  No  for  GA  and  N  =  qNG  for  DE. 


3.1.  Settings  for  the  Genetic  Algorithm 


GA  is  applied  to  E\  with  the  following  settings:  (q  =  64,/  =  0.8)  and  (q  =  128,/  =  0.4), 
both  with  Nc  =  500  (therefore  both  with  N  =  25600).  pc  was  varied  in  between  0. 1  and  0.9  in 
steps  of  0.2,  whereas  pm  was  varied  in  between  0.025  and  0.25  in  steps  of  0.025.  The  number 
of  independent  runs  Nr  equals  50  for  each  combination  of  settings.  The  estimates  of  the 
success  probabilities  for  the  considered  settings  arc  presented  in  the  figures  1  and  2. 
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Figure  1:  GA  with  q  =  64,  /  =  0.8;  statistical  error  on  ps  is  about  0.06 


Figure  2:  GA  with  q  =  128,  /  =  0.4;  statistical  error  on  ps  is  about  0.06 


It  is  observed  that  the  attainable  level  of  ps  hardly  depends  on  pc,  while  the  choice  of  pm  is 
of  major  importance.  There  is  a  rapid  fall  off  of  ps  versus  pm  for  pm  >  0.15.  The  optimum 
setting  of  pm  is  slightly  larger  than  0.1.  The  setting  q  =  128  and /  =  0.4  (also  with  JV=25600) 


leads  to  higher  levels  of  ps,  with  the  possibility  of  nearly  100  %  success.  Notice  the  need  of 
taking  the  (statistical)  error  on  ps  into  account!  Finally,  Gerstoft’s  default  setting  [1]  (q  =  64, 
Nc  =  500,/r  =  0.5,  pc  =  0.8,  pm  =  0.05)  is  not  robust,  leading  to  ps  =  0.69  (with  N  only  1 6000). 


3.2.  Settings  for  Differential  Evolution 

Applying  DE  to  E\,  two  combinations  of  q  and  Nc  are  considered:  (q  =  16,  Nc  =  150)  and 
(q  =  32,  Nc  =  75),  both  with  N  =  2400  (instead  of  the  25600  for  GA).  To  provide  more 
accurate  estimates  for  ps,  the  number  of  independent  runs  per  setting  is  taken  Nr  =  500.  ps  for 
different  combinations  of  F  and  pc  is  given  in  figures  3  and  4. 

It  is  seen  that  the  smaller  number  of  generations  Nc  =  75  requires  a  lower  setting  of  the 
contraction  factor  F.  DE  seems  to  be  more  sensitive  to  the  choice  of  F  than  to  that  of  the 
crossover  probability  pc. 


P«  pc 

Figure  3:  DE  with  q  =  16,  Ng  =  150  Figure  4:  DE  with  q  =  32,  Nc  =  75 
90%  ps  contour  97  %  ps  contour 


3.3.  Comparing  the  performance  of  GA  and  DE 

Table  1  gives  the  Gerstoft  GA  setting  [1]  and  the  optimum  settings  found  for  GA  and  DE. 


H 

Nc 

Pc 

Pm 

fr 

F 

N 

Ps 

GA  Gerstoft 

64 

500 

0.8 

0.05 

0.5 

16000 

0.69  ±  0.05 

GA  best 

128 

500 

0.9 

0.1 

0.4 

25600 

0.98  ±  0.02 

DE 

16 

150 

0.85 

0.8 

2400 

0.91  ±0.013 

DE 

32 

75 

0.75 

0.5 

2400 

0.98  ±  0.007 

DE 

32 

150 

0.8 

0.6 

4800 

0.99  ±  0.005 

Table  1:  GA  and  DE  best  settings  and  corresponding  performances. 

With  comparable  robustness  (ps  =  0.98)  DE  requires  only  2400  function  evaluations 
instead  of  the  25600  of  GA,  an  order  of  magnitude  improved  efficiency. 

How  to  explain  this  better  performance?  We  suspect  that  the  partner  creation  process  of 
DE  is  a  better  way  to  escape  from  a  local  minimum.  Therefore,  the  contraction  factor  F  is  an 
important  setting  parameter  of  DE.  At  the  same  time  this  partner  creation  process  offers  the 


possibility  to  better  explore  the  parameter  search  space,  provided  the  mutual  distances  of  the 
population  members  are  sufficiently  large. 


4.  SUMMARY  AND  CONCLUSIONS 

In  this  paper  a  relatively  new  global  optimisation  method  called  differential  evolution 
(DE)  has  been  tested.  This  method  can  be  applied  for  solving  inversion  problems  in 
underwater  acoustics.  Its  performance  in  terms  of  efficiency  and  robustness  was  compared 
with  that  of  a  genetic  algorithm  (GA).  Efficiency  is  measured  by  the  number  of  function 
evaluations.  Robustness  is  defined  as  the  probability  of  success.  The  fraction  of  runs  that 
converge  to  the  global  minimum,  out  of  a  number  of  independent  repetitions  of  the  algorithm 
for  a  given  setting  (typically  50),  estimates  this  probability. 

The  optimisation  problem  comprises  the  minimisation  of  an  energy  function  depending  on 
6  parameters,  featuring  various  local  minima,  parameter  coupling  and  varying  sensitivity  to 
the  parameters.  Minimising  this  function  with  GA  requires  several  tens  of  thousands  of 
function  evaluations.  The  performance  of  GA  hardly  depends  on  crossover  probability, 
whereas  it  depends  critically  on  mutation  probability.  Using  25600  function  evaluations 
almost  1 00  %  success  probability  was  realised,  provided  a  mutation  probability  of  1 0  %  is 
taken.  With  DE  the  same  success  rate  is  obtained  with  10  times  less  function  evaluations.  The 
performance  of  DE  depends  mainly  on  the  multiplication  factor  and  less  on  crossover 
probability. 

The  following  question  remains.  Do  the  best  settings  found  for  the  optimisation  algorithms 
depend  on  the  optimisation  problem  at  hand?  To  this  end  we  introduce  another  optimisation 
problem,  i.e.,  geo-acoustic  parameter  estimation  using  matched  field  inversion  [Snellen, 
2004]. 
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Abstract  As  the  vertical  directivity  of  ambient  noise  is  dependent  on  the  reflective  prop¬ 
erties  of  the  seabed,  geo-acoustic  parameters  can  he  inferred  by  inversion  of 
the  noise.  This  is  clearly  demonstrated  in  a  paper  by  (Harrison  and  Simons. 
2002).  In  this  work  the  seabed  reflection  loss  is  directly  found  by  comparing 
the  upward-  and  downward-going  noise.  These  are  determined  by  beamform¬ 
ing  the  ambient  noise  field  received  on  a  vertical  array  of  hydrophones.  Geo¬ 
acoustic  parameters  were  derived  with  a  hand-searching  technique  that  relies  on 
inspection  of  the  experimental  (angle-  and  frequency  dependent)  reflection  loss 
for  critical  angle,  fringe  spacing  and  modulation  depth.  In  the  current  paper 
by  applying  a  global  optimisation  technique  as  an  automatic  search  method  im- 
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proved  estimates  for  the  geo-acoustic  parameters  were  obtained.  For  the  global 
search  method  a  relatively  new  technique,  called  "differential  evolution",  is  ap¬ 
plied  (Storn  and  Price,  1995).  Compared  to  a  genetic  algorithm  a  ten  limes 
higher  efficiency  is  obtained,  since  differential  evolution  finds  the  optimal  pa¬ 
rameter  setting  with  only  one  tenth  of  the  number  of  forward  model  evaluations. 
For  these  inversions  two  different  models  for  the  seabed  reflection  loss  were  em¬ 
ployed.  First  we  assumed  the  seabed  to  consist  of  a  single  homogeneous  fluid 
sediment  layer  on  top  of  a  homogeneous  fluid  subbotlom.  As  the  fringe  pattern 
of  some  of  the  measured  reflection  loss  data  suggests  the  existence  of  at  least 
two  layers,  the  data  were  also  inverted  using  a  seabed  model  comprising  two 
llu id  sediment  layers  on  top  of  a  homogeneous  fluid  subbottom. 


Introduction 

As  the  coherence  and  vertical  directivity  of  ambient  noise  depends  on  the 
reflective  properties  of  the  seabed,  geo-acoustic  parameters  can  be  inferred  by 
inversion  of  the  noise.  This  is  clearly  demonstrated  by,  among  others,  (Car¬ 
bone  et  al.t  1998)  and  (Harrison  and  Simons,  2002). 

Carbone  measures  the  broadband  noise  coherence  with  a  vertically  separated 
pair  of  hydrophones,  which  was  compared  with  noise  model  calculations  search¬ 
ing  over  geo-acoustic  parameter  space. 

In  (Harrison  and  Simons,  2002)  an  approximation  of  the  plane  wave  reflec¬ 
tion  coefficient  of  the  seabed  (as  a  function  of  grazing  angle  and  frequency)  is 
directly  found  by  comparing  the  upward-  and  down  ward- going  noise.  These 
are  determined  by  beamforming  the  ambient  noise  field  received  on  a  vertical 
array  of  hydrophones.  Geo-acoustic  parameters  were  derived  with  a  hand¬ 
searching  technique  that  relies  on  inspection  of  the  experimental  (angle-  and 
frequency  dependent)  reflection  loss  for  critical  angle,  modulation  depth  and 
fringe  spacing.  From  these  phenomena  sound  speed,  impedance  mismatch 
and  layer  thickness  can  already  be  determined  fairly  well.  The  potential  of 
this  method  was  demonstrated  with  experimental  data  from  several  sites  in  the 
Mediterranean,  Pros  and  cons  of  the  method  are  extensively  discussed.  Im¬ 
portant  strengths  of  the  method  are  that  no  noise  model  is  required  and  that 
it  tolerates  an  unknown  noise  source  distribution  and  noise  from  distant  ship¬ 
ping.  The  paper  of  (Harrison  and  Simons,  2002)  does  not  go  further  than  hand 
searches  for  the  geo-acoustic  parameters,  since  the  main  point  was  to  demon¬ 
strate  merely  that  a  solution  is  possible.  It  was  concluded  that  by  applying  a 
global  optimisation  technique  as  an  automatic  search  method,  improved  esti¬ 
mates  for  the  geo-acoustic  parameters  should  be  obtained.  It  should  however 
be  emphasised  that  results  from  this  inversion  are  not  necessarily  unique, 
in  this  paper  we  compare  the  measured  reflection  coefficient  with  calculations 
using  models  for  the  plane  wave  reflection  coefficient,  searching  over  geo¬ 
acoustic  parameter  space.  For  these  inversions  we  first  assume  the  seabed  to 
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consist  of  a  single  homogeneous  fluid  sediment  layer  on  top  of  a  homogeneous 
fluid  subbottom.  As  the  fringe  pattern  of  some  of  the  measured  reflection  loss 
data  suggests  the  existence  of  at  least  two  layers,  the  data  are  also  inverted 
using  a  seabed  model  comprising  two  fluid  sediment  layers  on  top  of  a  homo¬ 
geneous  fluid  subbottom.  These  models  are  described  in  section  1.  Also,  the 
smearing  effect  of  the  beamforming  is  discussed  in  that  section. 

For  the  global  search  method  a  relatively  new  technique,  called  "differential 
evolution"  (Stom  and  Price,  1995),  is  applied.  A  brief  description  of  differen¬ 
tial  evolution,  and  a  comparison  with  a  genetic  algorithm,  is  given  in  section  2. 
The  ambient  noise  experiments  and  the  data  selected  for  this  study  are  de¬ 
scribed  in  section  3.  The  inversion  results  are  given  and  discussed  in  section  4. 
Section  5  presents  a  summary  and  some  conclusions. 

1.  Ambient  noise  inversion  -  applied  theory 

A  beam -smeared  approximation  of  the  bottom  reflection  coefficient  is  ob¬ 
tained  by  calculating  the  ratio  of  the  downward  and  upward  array  response. 
Theoretically,  this  is  given  by  (Harrison  and  Simons,  2001) 


Aj-i %)  f-j%  /)  2 9M 

A{6s)  /-f/2  T=mnm  e”  /> sin  2m 


with 

A  +.  upward  array  response 

A  _  downward  array  respoase 

9  grazing  angle  of  incidence  at  the  seabed 

0*  steering  angle  of  the  beam  former 

/  frequency 

R  true  bottom  reflection  coefficient 

B  beam  pattern  of  the  vertical  array 


The  real  bottom  loss  BL  is  given  by  BL  =  —20  log]0 1  R\,  whereas  the  beam- 
smeared  bottom  loss  is  given  by  BL  =  — 201ogm  |i?|.  The  effect  of  beamform¬ 
ing  is  to  smear  out  features  in  BL,  especially  at  the  lower  frequencies. 

We  search  with  an  optimisation  method  for  the  geo-acoustic  parameters  by 
matching  measured  and  modelled  bottom  loss  BL  using  the  following  energy 
function 


E(m)  =  ~"£(BLcgl(9J,m)  -  BLcxp(6J))2  (1.2) 


with 
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Figure  1.1a.  The  single-sediment 
geo-acoustic  bottom  model. 


Figure  1.1b.  The  two-sediment 
geo-acoustic  bottom  model. 


m  the  vector  containing  the  unknown  geo-acoustic  parameters 

C  a  normalisation  constant 

BLcai  calculated  or  modelled  bottom  loss  in  dB  (including  smearing  effect) 
BLtXp  experimental  or  measured  bottom  loss  in  dB 


Two  models  for  the  seabed  have  been  employed:  a  seabed  consisting  of 
a  single  fluid  sediment  layer  on  top  of  a  homogeneous  fluid  subbottom  (see 
figure  1.1a)  and  a  seabed  consisting  of  two  fluid  homogeneous  sediment  layers 
on  top  of  a  fluid  subbottom  (see  figure  1  .lb). 

Each  sediment  and  subbottom  is  assumed  to  have  a  constant  sound  speed, 
density  and  attenuation.  To  limit  the  size  of  the  parameter  search  space,  we  do 
not  invert  for  the  densities  p2,  ps  and  Pa.  Density  is  coupled  to  sound  speed 
according  to  one  of  the  empirical  relations  given  by  (Hamilton  and  Bachman, 
1982).  Further,  we  assume  <24  =  03  =  a-i  =  a,  because  the  reflection  co¬ 
efficient  is  less  sensitive  to  attenuation  in  the  deeper  layers.  For  these  layered 
structures  analytical  solutions  for  the  reflection  coefficient  are  given  by  (Jensen 
etal.,  1994). 

For  the  two-sediment  layered  structure  the  reflection  coefficient  R  is  given 
by 


R12  +  R-24  exp(2t<fo) 
1  +  RV2R24  exp(2i02) 


with 


R‘23  +  Ru  exp(2i<fe) 
1  +  R-iaRm  exp  (2 ) 


(1.4) 


where  Rjk  is  the  (Rayleigh)  reflection  coefficient  for  the  interface  between 
medium  j  and  k  (medium  1  =  water,  medium  2  =  sediment  no.  I  ,  etc.,  see 
figure  1): 
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PkCk  _  PjCj 
D  _  sine*  sinfl, 

KJk  ~  ImZ  |  E  5~ 

sin  Ok  '  sin  B0 


(1.5) 


for  k  =  j  +  1;  j  =  1,2, 3.Here  8  =  6\  and  the  vertical  plane  phase  delays 
are  given  by  =  /r2/i2sm02  and  fo  =  k^h^sindi. 

Further,  Snell’s  law  applies:  k^cosOi  =  k^cosOs  =  k^cosO^  =  kicosO i 
where  the  wavenumber  in  medium  j  is  given  by  kj  = 

For  the  single-sediment  situation  the  reflection  coefficient  is  given  by 


Rn  +  fl23exp(2i02) 

1  +  f?i2f?23exp(2i02) 

In  these  equations  the  sound  speed  Cj  is  complex  when  the  attenuation  aj  ^ 
0.  The  imaginary  part  of  Cj  is  then  given  by 

- CjOj -  =  _CjOj_ 

40nlogU)(e)  54.575  '  ; 


2.  Description  of  differential  evolution 

In  underwater  acoustics  genetic  algorithms  (GA)  were  introduced  by  (Ger- 
stoft,  1994).  Its  principle  of  operation  is  summarised  as  follows.  To  optimise 
an  energy  function  E .  GA  improves  populations  of  q  parameter  combinations 
during  successive  generations.  A  GA  optimisation  run  starts  by  creating  at  ran¬ 
dom  an  initial  population  of  size  q.  Each  member  represents  a  certain  parame¬ 
ter  value  combination  m  and  its  energy  E  is  calculated.  Because  a  population 
member  with  a  lower  E  has  a  higher  fitness,  it  attains  a  higher  probability  to 
be  selected  in  one  (or  more)  of  the  pairs  of  parents  to  be  formed.  This  results 
in  q/2  pairs  of  parents  with  a  larger  proportion  of  fit  members.  From  each  pair 
of  parents  two  children  are  created  by  a  random  exchange  of  parental  values, 
steered  by  the  crossover  probability  pc,  followed  by  random  changes  of  the  Np 
parameter  values,  with  mutation  probability  pm.  A  set  of  q  children  vectors 
results.  The  next  generation  is  established  by  taking  at  random  frq  members 
of  these  children  and  selecting  the  (1  —  fr)q  most  fit  members  of  the  current 
population.  This  process  is  repeated  over  several  hundreds  of  generations. 

A  GA  has  the  following  setting  parameters:  population  size  q ,  crossover  proba¬ 
bility  pr,  mutation  probability  pm,  reproduction  size  fr  and  the  number  of  gen¬ 
erations  Nc ■  These  have  to  be  selected  carefully  and  can  be  problem-specific, 
i.e..  they  depend  on  the  optimisation  problem  at  hand  (see  e.g.  (van  Moll  and 
Simons,  2004)). 

Differential  evolution  (DE).  just  like  GA.  starts  with  an  initial  population  of  q 
randomly  chosen  parameter  value  combinations  m.  These  m’s  are  improved 
during  successive  generations  of  constant  size  q ,  in  the  sense  that  a  descendant 
replaces  an  ni.  becoming  its  successor,  if  it  has  a  lower  energy.  The  distinctive 
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feature  of  DE  is  the  way  in  which  these  descendants  are  created.  In  (Storn  and 
Price,  1995)  various  ways  to  generate  new  m's  are  described.  Here  only  the 
following  procedure  is  considered.  At  the  start  of  generation  k  the  parameter 
vectors  m*  7  are  given  and  for  each  of  them  a  descendant  is  created, 

being  a  potential  successor.  To  create  this  descendant  x  of  first  a  part¬ 
ner  is  constructed  according  to 

PM  =  mkiil  +  F(mk,i2  -  nujj)  (1-8) 

with  three  different  vectors  m^,,  m*,j3  and  m^i3  chosen  at  random  from 
the  population  and  F  a  scalar  multiplication  factor  between  0  and  1 .  The  de¬ 
scendant  dkii  of  mu  results  from  applying  crossover  to  nn,,  and  p*,j  with 
crossover  probability  pc.  A  higher  value  of  p,-  leads  (on  the  average)  to  more 
dimensions  of  p k,i  to  be  copied  into  dfc,,,  while  the  values  for  the  remaining 
dimensions  are  taken  from  mkti3.  d*,*  only  replaces  becoming  its  suc¬ 
cessor.  if  its  energy  is  lower.  The  setting  parameters  of  DE  are  population  size 
q.  multiplication  factor  F,  crossover  probability  pc  and  the  number  of  genera¬ 
tions  Ng. 

Comparing  GA  and  DE,  both  show  improvement  of  a  fixed  sized  population  of 
parameter  combinations.  However,  in  each  generation  of  GA  a  fixed  fraction 
of  the  population  is  replaced  by  children,  whether  they  are  better  or  not.  In  DE 
each  member  of  the  population  will  be  replaced  by  its  descendant  only  if  this 
descendant  has  lower  energy.  On  the  other  hand,  the  fitness  of  a  member  in¬ 
creases  the  probability  of  that  member  to  play  a  role  in  the  creation  of  children 
in  GA,  where  in  DE  this  fitness  does  not  influence  the  probability  of  a  mem¬ 
ber  to  be  selected.  Another  difference  is  the  way  a  partner  p*,j  is  found  for 
irifej.  While  in  GA  the  partner  is  an  element  of  the  given  generation,  in  DE  it 
is  constructed  from  the  population  and  thus  not  an  element  of  this  generation. 
The  partner  creation  process  in  DE  is  steered  by  the  multiplication  factor  F  and 
replaces  the  mutation  mechanism  of  GA. 

In  (van  Moll  and  Simons,  2004)  and  (Snellen  et  al„  2004)  it  was  demonstrated 
that,  compared  to  a  GA,  a  ten  times  higher  efficiency  is  obtained  with  DE,  as 
DE  finds  the  global  optimum  with  only  one  tenth  of  the  number  of  forward 
model  evaluations. 

3.  The  selected  ambient  noise  data 

For  this  study  we  selected  the  North  Elba  case  (silt  bottom,  water  depth  1 2 1 
m)  from  the  five  sites  in  the  Mediterranean  where  ambient  noise  measurements 
were  carried  out.  The  motivation  for  selecting  the  North  Elba  data  is  the  irreg¬ 
ularity  of  the  interfering  fringes  in  the  observed  bottom  loss,  which  might  be 
evidence  for  the  presence  of  two  sediment  layers.  This  is  considered  to  be  a 
challenging  situation  for  inversion.  The  experiments  were  described  in  great 
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detail  by  (Harrison  and  Simons,  2002).  They  were  carried  out  with  a  moored 
64  element  vertical  array  with  an  acoustic  aperture  of  62  m  roughly  centred  in 
the  water  column.  Only  the  middle  32  regularly  spaced  hydrophones  (aperture 
16  m,  element  spacing  0.5  m)  were  used  for  the  processing,  which  consisted 
of  beamforming  of  manually  selected  portions  of  noise  data.  Each  portion 
comprises  a  few  seconds  of  ambient  noise  data.  The  resulting  array  responses 
versus  angle  and  frequency  were  averaged  over  many  portions  of  data.  The 
up-to-down  ratio  of  the  array  response  is  thus  based  on  a  few  minutes  of  data 
in  total. 


4.  Inversion  results 

4.1  Single  sediment  layer  inversions  -  the  North  Elba  data 

We  discuss  the  inversion  results  of  the  bottom  loss  data  obtained  from  the 
North  Elba  ambient  noise  data  (see  figure  1,5a)  assuming  a  single  sediment 
model  for  the  seabed  (see  figure  1.1a).  As  explained  in  section  1  we  do  not 
invert  for  density,  hence  the  number  of  parameters  sought  amounts  to  4.  The 
search  bounds  of  the  parameters  to  be  inverted  are  given  in  table  1.1.  A  con¬ 
straint  in  the  inversion  algorithm  is  C3  >  02  (cj  =  1514  m/s). 


Table  1.1.  Parameters  to  be  inverted  for:  single -sediment  inversions  of  the  North  Elba  data. 


Parameter 

Search  range 

Inversion  result 

C2  [m/sl 

1400-1700 

1520.4  ±  0.6 

pi  [g/cma  ] 

RA 

1.562  ±0.002 

hi  [m| 

0.1-5 

3.40  ±  0.03 

a  [dB/Aj 

0-2 

0.00031  ±  0,00004 

c3  [m/s| 

1500-1900 

1538.4  ±  1.5 

pa  Ig/cm3] 

RA* 

1.020  ±  0.004 

DE  was  run  57  times.  In  the  following  we  perform  a  statistical  analysis  of 
the  final  parameter  estimates  of  52  successful  runs.  Parameter  correlation  is 
illustrated  in  figure  1 .2  showing  the  magnitude  of  the  linear  correlation  coeffi¬ 
cient  calculated  for  all  combinations  of  the  parameter  estimates.  The  statistical 
significance  of  each  of  the  correlation  coefficients  was  also  determined:  the  red 
stars  in  the  figure  indicate  for  which  parameter  combination  the  confidence  of 
the  observed  correlation  coefficient  exceeds  95%. 

In  figure  1.3  parameter  estimates  are  plotted  against  each  other  for  combi¬ 
nations  with  strong  coupling.  There  is  no  correlation  between  the  sediment 
geo-acoustic  parameters  (c^,  C3  and  a)  and  the  sediment  thickness  hi. 

Histograms  of  the  parameter  estimates  are  given  in  figure  1.4.  The  corre¬ 
sponding  mean  and  standard  deviations  for  each  parameter  are  given  in  table 
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Figure  L2.  Single- sediment  inversions  of  the  North  Elba  data:  magnitude  of  the  correlation 
coefficient  for  all  parameter  combinations,  A  red  star  is  plotted  if  the  correlation  is  statistically 
significant. 


(IiWb) 


Figure  1 J*  Single-sediment  inversions  of  the  North  Elba  data:  The  52  parameter  estimates 
plotted  against  each  other  for  those  combinations  with  strong  coupling. 


1.1.  It  can  be  concluded  that  all  parameters  are  accurately  determined.  The  pa¬ 
rameter  values  obtained  are  further  assessed  in  section  4.2.2  where  we  discuss 
inversions  of  the  same  data  with  a  two-sediment  layer  model. 

Using  equation  1.1  we  have  calculated  the  beam-smeared  bottom  loss  BL 
for  the  mean  values  of  the  inversion  results  as  given  in  table  1.1,  the  result 
of  which  is  given  in  figure  1.5d.  The  deviation  from  the  measured  (beam- 
smeared)  bottom  loss,  as  presented  in  figure  1.5a,  is  given  in  figure  1.5e.  The 
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Figure  L4.  Single- sediment  inversions  of  the  North  Elba  data:  histograms  of  the  52  estimates 
for  each  parameter 


same  is  done  for  the  parameters  reported  in  (Harrison  and  Simons,  2002), 
which  were  obtained  using  the  hand-searching  technique  (see  figures  1 .5b  and 
1.5c).  We  observe  that  a  considerable  improvement  is  obtained  with  an  auto¬ 
mated  search  technique 

4.2  Inversions  with  two  sediments 

4.2.1  Synthetic  case.  First,  the  performance  of  the  inversion  with  two 
sediments  is  assessed  using  a  synthetic  reflection  coefficient  created  for  the 
parameter  values  given  in  table  1,2,  Again,  we  do  not  invert  for  density  (see 
section  1),  hence  the  number  of  parameters  sought  amounts  to  6.  Constraints  in 
the  inversion  algorithm  are:  c4  >  max(c2 ,  C3)  and  h 2  4-  h$  <  5  m  (ci  =  1514 
m/s). 

The  inversion  was  run  171  times.  On  the  basis  of  the  energy  function  values 
obtained  for  these  runs,  two  solutions  to  the  inversion  problem  can  be  distin¬ 
guished,  see  figure  1.6. 

First,  we  select  the  best  solutions  for  which  E  <  0.01  (111  out  of  171).  For 
these  solutions  parameter  correlation  is  illustrated  in  figure  1.7,  which  shows 
the  magnitude  of  the  linear  correlation  coefficient  calculated  for  all  combina¬ 
tions  of  parameter  estimates.  Again,  w^e  indicate  for  which  parameter  combi¬ 
nation  the  confidence  of  the  observed  correlation  coefficient  exceeds  95  %,  As 
in  the  case  of  the  single  sediment  inversions  of  measured  data  (see  previous 
paragraph),  the  geo-acoustic  parameters  of  the  two  sediments  (02,03,04  and 
a)  are  strongly  coupled.  Also,  the  thicknesses  of  the  sediments,  hi  and  h 3*  are 
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Figure  L5a. 

Measured  bottom  loss 
BL 


Figure  JJb. 

Modelled  BL  using  geo 
acoustic  parameters  from 
(Harrison  and  Simons.  2002) 


Figure  L5d 
Modelled  BL  using 
DE- inverted  geo -acoustic 
parameters  (single  sediment) 


Figure  Uf. 
Modelled  BL  using 
DE- inverted  geo- 
acoustic  parameters 
{two  sediments) 


Figure  1.5c * 

Deviation  of  modelled 
BL  (figure  1,5b) 
from  measured  BL 


Figure  L5e. 

[deviation  of  modelled 
BL  (figure  L5d) 
from  measured  BL 


} 


i,.*'  j 

Hi  40  B  #0 

(FUinfl  #r.yb*  |'*g| 


Figure  1.5g. 

Deviation  of  modelled 
BL  (figure  L5f) 
from  measured 
BL 
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Table  1.2.  Parameters  to  be  inverted  for:  two-sediment  inversions  of  the  synthetic  data. 


Parameter 

Search  range 

True  value 

Inversion  result 

Best  solution  E  <  0.01 

Inversion  result 
Solutions  with  E  &  0.07 

C2  [m/s] 

1400-1700 

1540 

1540.1  ±0.3 

1539.1  ±0,2 

pi  [g/cm:,l 

N.A 

±624 

1.625  ±0.001 

±622  ±0.001 

hi  [m] 

0,1-5 

2.0 

2.000  ±0.007 

3.006  ±  0.002 

cs  (m/s] 

1400- 1700 

1530 

1530.1  ±0,3 

1575.1  ±0.3 

p»  Ig/cnr'l 

N.A 

1.594 

1.595  ±0.001 

1,714  ±0.001 

h3  Iml 

0.1-5 

1.0 

1.001  ±  0,007 

1.3  ±0.4 

('4  I  m/s} 

1500-1900 

1570 

1570,1  ±0.4 

1575.1  ±0,2 

Pi  Ig/cm3! 

N.A 

1.702 

1,702  ±0.001 

±714  ±0.001 

a  IdBM| 

0-2 

0.05 

0.051  ±0,001 

0.055  ±0.001 

Figure  1.6.  Two-sediment  inversions  of  the  synthetic  data:  histogram  of  the  final  energy  func¬ 
tion  values  obtained  for  the  17 1  successful  inversions. 


strongly  coupled,  but  there  is  virtually  no  correlation  between  the  geo-acoustic 
parameters  and  the  thicknesses. 

Histograms  of  these  parameter  estimates  are  given  in  figure  1 .8.  The  cor¬ 
responding  mean  and  standard  deviation  for  each  parameter  are  given  in  table 
1 .2.  It  can  be  concluded  that  all  parameters  are  accurately  determined  and  that 
their  mean  values  are  in  excellent  agreement  with  the  true  values. 

Next,  we  select  the  somewhat  worse  solutions  where  E  is  around  0.07  (see 
the  second  peak  in  the  histogram  for  E,  figure  1 .6).  Histograms  of  these  param¬ 
eter  estimates  are  given  in  figure  1.9.  The  corresponding  mean  and  standard 
deviations  for  each  parameter  are  given  in  table  1.2.  It  is  observed  that  a  sin¬ 
gle  layer  solution  is  found  now.  i.e.,  a  top  sediment  layer  with  a  sound  speed 
nearly  equal  to  that  of  the  first  layer,  with  a  thickness  equal  to  the  sum  of  the 
two  sediments  (being  3  m).  Further,  C3  becomes  equal  to  C4  and.  consequently. 
/13  becomes  undetermined.  Still,  this  somewhat  worse  solution  is  physically 
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Figure  1.7.  Two-sediment  inversions  of  the  synthetic  data:  magnitude  of  the  correlation  co¬ 
efficient  for  all  parameter  combinations,  A  red  star  is  plotted  if  the  correlation  is  statistically 
significant. 


Figure  1.8.  Two-sediment  inversions  of  the  synthetic  data:  histograms  of  the  best  estimates 
for  each  parameter  ( E  <  GT)I). 
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Figure  1.9.  Two-sedimem  inversions  of  ihe  synthetic  data:  histograms  of  the  somewhat  worse 
estimates  ( E  «  0.07). 


plausible.  The  somewhat  higher  sound  speed  found  in  the  lower  !  m  of  the  3 
m  sediment  (1540  m/s  instead  of  the  true  1530  m/s)  and  also  below  the  sedi¬ 
ment  (1575  m/s  instead  of  the  true  1570  m/s)  is  compensated  by  a  somewhat 
higher  attenuation  (0.055  dB/A  instead  of  the  true  0.05  dB/A ). 

The  good  performance  of  this  synthetic  inversion  problem  is  promising  for 
inversion  of  measured  data  using  a  two-sediment  geo-acoustic  model.  How¬ 
ever.  the  two  solutions  discussed  also  illustrate  the  non-uniqueness  of  the  solu¬ 
tion  of  geo-acoustic  inverse  problems.  For  this  synthetic  case  the  two  solutions 
can  clearly  be  distinguished  by  their  energy.  For  the  inversion  of  experimental 
data  it  is  expected  that  it  can  be  less  obvious  which  solution  to  choose. 

4.2.2  The  North  Elba  data.  We  now  discuss  the  inversion  results  of 
the  bottom  loss  data  obtained  from  the  North  Elba  ambient  noise  data  (see  fig¬ 
ure  1.5a)  assuming  a  two-sediment  model  for  the  seabed  (see  figure  I.  lb).  The 
inversion  algorithm  is  the  same  as  that  for  the  synthetic  inversions  described 
in  the  previous  paragraph.  The  inversion  was  run  200  times,  of  which  1 14 
successful  runs  were  selected.  Histograms  of  the  corresponding  parameter  es¬ 
timates  are  given  in  figure  1.10.  It  seems  that  sediment  thicknesses  are  not  well 
resolved.  However,  for  the  total  sediment  thickness,  h>2  +  h$  ,  we  do  obtain  a 
peaked  distribution.  Also,  as  in  the  synthetic  case,  h 2  and  ha  are  correlated. 

Figure  1.11  presents  a  plot  of  the  estimates  for  h 2  plotted  against  those  for 
ha-  It  seems  obvious  to  make  a  selection  based  on  total  sediment  thickness. 
We  chose  3.3  <  h-2  +  hs  <  3.5  m. 

Histograms  of  the  selected  solutions  are  given  in  figure  1.12.  (Of  course, 
the  distribution  for  h-2  +  h$  is  very  peaked  now,  but  /12  and  ha  separately  re- 
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Figure  1.  IQ.  Two-sediment  inversions  of  the  North  Elba  data:  histograms  of  the  U4  estimates 
for  each  parameter  Note  that  we  included  the  histogram  for  total  sediment  thickness. 


Figure  JJ1.  Two-sediment  inversions  of  the  North  Elba  data:  114  estimates  for  h2  plotted 
against  those  for  ha*  The  selected  solutions  are  indicated  by  the  red  symbols. 
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Figure  1.12.  Two-sediment  inversions  of  the  North  Klba  data:  histograms  of  the  parameter 
estimates,  whieh  were  selected  on  total  sediment  thickness. 


main  unresolved).  The  corresponding  mean  and  standard  deviations  for  each 
parameter  are  given  in  table  1.3.  The  inversion  result  reported  in  (Harrison  and 
Simons.  2002)  is  also  included  for  comparison. 


Table  1.3.  Parameters  to  be  inverted  for:  two-sediment  inversions  of  the  North  Ivlba  data. 


Parameter 

Search  range 

Inversion  result 
Harrison  and  Simons 

Inversion  result 
This  work 

t'2  I  m/s] 

1400-1700 

1600 

1523  ±  4 

P2  [g/cm3] 

N.A 

1.8 

1.57  ±0.01 

h2  |m] 

0.1-5 

0.7 

1.57  ±  0.65 

Cs  [m/s] 

1400-1700 

1530 

1523  ±  7 

p:\  |g/cm3] 

N.A 

1.75 

1.57  ±0.03 

U:\  [m] 

0.1-5 

2.1 

1.82  ±0.65 

ca  1  m/s] 

1500-1900 

1600 

1546  ±  8 

(M  [g/cm3] 

N.A 

1.8 

1.64  ±0.02 

a  IdB/A] 

0-2 

0.141 

0.0005  ±  0.0003 

h2  +  /i3  [m] 

N.A. 

N.A. 

3.38  ±  0.03 

For  completeness,  parameter  coupling  of  the  selected  solutions  is  shown  in 
figure  1 . 1 3.  a  similar  behaviour  as  that  of  the  synthetic  inversion  problem  (fig¬ 
ure  1 .7).  Basically,  a  single  sediment  solution  is  obtained,  because  the  inver¬ 
sion  result  for  c-2  equals  that  of  C3  and  a  peaked  distribution  is  only  obtained  for 
h'2  -F  /13  and  not  for  hi  or  /13  separately.  Further,  the  solution  agrees  well  with 
the  single-sediment  solution  given  in  table  1.1.  However,  it  does  not  agree  well 
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Figure  1.13.  Two-sediment  inversions  of  the  North  Elba  data:  magnitude  of  the  correlation 
coefficient  for  the  selected  solutions.  A  red  star  is  plotted  if  the  correlation  is  statistically  sig¬ 
nificant. 


with  the  solution  of  (Harrison  and  Simons,  2002),  albeit  that  both  solutions  are 
representative  of  a  silt  bottom. 

Finally,  we  have  calculated  the  beam-smeared  bottom  loss  BL  for  the  mean 
values  of  the  inversion  results  as  given  in  table  1.3,  the  result  of  which  is  given 
in  figure  I.5f.  The  deviation  from  the  measured  (beam-smeared)  bottom  loss, 
as  presented  in  figure  1,5a,  is  given  in  figure  l,5g.  No  further  improvement 
compared  to  the  single-sediment  solution  (figures  1.5d  and  L5e)  is  obtained, 
wfhich  is  not  surprising,  as  the  single-sediment  and  two-sediment  inversion  re¬ 
sults  are  basically  the  same. 

5.  Summary  and  conclusions 

An  approximation  of  the  plane  wave  reflection  coefficient  of  the  seabed  (as 
a  function  of  grazing  angle  and  frequency)  can  be  found  by  comparing  the 
upward-  and  downward-going  noise.  In  (Harrison  and  Simons,  2002)  these 
were  determined  by  beamforming  the  ambient  noise  field  received  on  a  verti¬ 
cal  array  of  hydrophones.  Geo-acoustic  parameters  were  derived  with  a  hand- 
searching  technique  that  relies  on  inspection  of  the  experimental  reflection  loss 
for  critical  angle,  modulation  depth  and  fringe  spacing.  In  this  article  we  in¬ 
ferred  geo-acoustic  parameters  by  inversion  of  the  measured  reflection  coef¬ 
ficient  with  a  global  optimisation  technique  as  an  automatic  search  method. 
For  this  study  we  selected  the  North  Elba  data  (silt  bottom),  one  of  the  sites 
in  the  Mediterranean  where  ambient  noise  measurements  were  carried  out  (see 
(Harrison  and  Simons,  2002)).  The  measured  reflection  coefficient  was  com- 


Inversion  of  shallow  water  ambient  noise  data 
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pared  with  calculations  using  models  for  the  plane  wave  reflection  coefficient, 
searching  over  geo-acoustic  parameter  space.  For  these  inversions  we  assumed 
the  seabed  to  consist  of  a  single  homogeneous  fluid  sediment  layer  on  top  of 
a  homogeneous  fluid  subbottom.  The  data  were  also  inverted  using  a  seabed 
model  comprising  two  fluid  sediment  layers  on  top  of  a  homogeneous  fluid 
subbottom.  For  the  global  search  method  a  relatively  new  (but  fast)  technique, 
called  "differential  evolution”  (Stom  and  Price.  1995),  was  applied.  It  was 
found  that  all  parameters  could  be  accurately  determined.  With  the  two-layer 
geo-acoustic  model  basically  a  single-sediment  solution  is  obtained,  which  is 
in  excellent  agreement  with  the  solution  obtained  with  the  single-sediment 
model.  However,  they  do  not  agree  well  with  the  solution  reported  in  (Har¬ 
rison  and  Simons,  2002),  albeit  that  all  solutions  are  representative  of  a  silt 
bottom.  A  comparison  of  the  beam-smeared  bottom  loss,  calculated  for  the 
geo-acoustic  parameters  inverted,  with  the  measured  bottom  loss  shows  that 
a  considerable  improvement  can  be  obtained  with  an  automated  search  tech¬ 
nique. 

Notes 

1 .  Attenuation  of  firs*  sediment  layer 
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Abstract;  The  angle  and  frequency  dependent  reflective  properties  of  the  seabed  can  be 
determined  from  beam-steered  ambient  noise  measurements  on  a  vertical  array  of 
hydrophones.  From  the  up-to-down  ratio  the  beam-smeared  modulus  square  of  the  plane 
wave  reflection  coefficient  is  obtained.  Geo-acoustic  parameters  can  be  inverted  by  matching 
the  measured  and  modelled  reflection  loss  assuming  a  certain  number  of  sediment  layers.  The 
searching  for  the  geo-acoustic  parameters  is  done  with  fast  global  optimisation  techniques 
such  as  differential  evolution.  An  alternative  approach  is  to  Fourier  transform  the  reflection 
coefficient  producing  the  minimum  phase  impulse  response,  which  explicitly  shows  layer 
boundaries ,  Spectral  factorisation  is  applied  to  restore  the  lost  phase  of  the  reflection 
coefficient.  This  method  offers  the  poss  ib  ility  ’  of  determining  the  number  of  layers  and  their 
order.  Both  methods  were  applied  to  reflection  loss  data  extracted  from  ambient  noise 
measurements  on  a  drifting  vertical  array  in  the  shallow  waters  in  the  Mediterranean.  The 
results  of  both  methods  are  compared  showing  good  agreement  between  the  two  sets  of 
sediment  layering  estimates.  It  is  concluded  that  the  methods  offer  a  sub-bottom  profiling 
technique  without  making  additional  sound  where  the  array  is  sensing  local  seabed 
properties  with  a  footprint  of  order  the  water  depth , 

Keywords;  Geo-acoustic  inversion ,  ambient  noise ,  sub-bottom  profiling 


1.  INTRODUCTION 

In  a  series  of  recent  papers,  methods  for  deriving  geo-acoustic  parameters  from  measured 
ambient  noise  fields  are  presented,  Carbone  ([!])  measures  the  broadband  noise  coherence 


with  a  vertically  separated  pair  of  hydrophones,  w'hich  is  compared  with  noise  model 
calculations  searching  over  geo-acoustic  parameter  space.  In  [2]  an  approximation  of  the 
plane  wave  reflection  coefficient  of  the  seabed  (as  a  function  of  grazing  angle  and  frequency) 
is  determined  directly  from  the  beam- formed  ambient  noise  field  received  on  a  vertical  array 
of  hydrophones  by  comparing  the  upward*  and  downw  ard-going  noise.  Subsequently,  model 
calculations  for  the  reflection  coefficient  and  an  efficient  global  search  technique  are 
employed  for  inferring  the  geo-acoustic  parameters  ([3]).  In  [4]  the  seafloor  layering  is 
obtained  in  a  direct  manner  by  Fourier  transforming  the  reflection  coefficient  as  a  function  of 
frequency  for  each  angle.  Spectral  factorization  is  applied  for  restoring  the  phase  of  the 
reflection  coefficient. 

In  this  paper  we  apply  the  latter  two  methods,  both  employing  the  measured  reflection 
coefficient,  for  sub-bottom  profiling.  To  demonstrate  the  methods'  capabilities  they  are 
applied  to  real  data  consisting  of  a  series  of  ambient  noise  data  acquired  writh  a  drifting 
vertical  hydrophone  array. 

In  section  2  the  method  of  inverting  the  reflection  coefficient  data  is  presented.  For  the 
global  search  the  technique  of  “differential  evolution”  ([5])  is  applied.  A  brief  description  of 
differential  evolution  is  also  included  in  section  2.  In  section  3  the  method  of  Fourier 
transforming  the  reflection  coefficient  is  described.  The  ambient  noise  experiments  and  the 
data  selected  for  this  study  are  described  in  section  4.  The  sub-bottom  profiling  results  as 
obtained  by  both  methods  arc  presented,  compared  and  discussed  in  section  5.  Section  6  gives 
a  summary  and  some  conclusions. 


2*  SLB  BOTTOM  PROFILING  THROUGH  INVERSION  OF  SHALLOW  WATER 
AMBIENT  NOISE 


A  beam-smeared  approximation  R  of  the  bottom  reflection  coefficient  R  is  obtained  by 
determining  the  ratio  of  the  downward  and  upward  array  response  of  the  vertical  array. 
Theoretically,  this  is  given  by  ([3]) 
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with  and  A_  the  upward  and  downward  array  response,  respectively.  6  is  the  grazing 
angle  at  the  seabed,  ft  the  steering  angle  of  the  be  am -former,  /  the  frequency,  R  the  true 
bottom  reflection  coefficient  and  b  the  beam  pattern  of  the  vertical  array. 

The  real  bottom  loss  B  is  given  by 
B  -  -20lCilog  |  R  |  and  B  =  -20lolog  |  R  |  is 
the  beam-smeared  bottom  loss.  The  effect  of 
beam-forming  is  to  smear  out  features  in  B , 
especially  at  the  lower  frequencies.  The 
seabed  is  assumed  to  consist  of  a  single  fluid 
sediment  layer  on  top  of  a  homogeneous 
fluid  sub-bottom  (see  Fig,  1),  The  sediment 
and  sub-bottom  are  each  assumed  to  have  a 
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Fig.  J:  The  single-sediment  geo-acoustic 
bottom  model . 


constant  sound  speed,  density  and  attenuation*  Further,  the  attenuation  in  the  sediment  and 
sub-bottom  arc  assumed  equal,  i.c*,  ai  = 

Analytical  solutions  for  the  reflection  coefficient  as  given  by  [6]  arc  used  for  the  forward 
calculations.  We  search  with  an  optimisation  method  for  the  six  unknown  geo-acoustic 
parameters,  m  =  ( a ,  pi,  h%,  p 3,  a),  by  matching  measured  and  modelled  bottom  loss  B 
using  the  following  energy  function 


(2) 


with  m  the  vector  containing  the  unknown  geo-acoustic  parameters,  C  a  normalisation 
constant,  5fa/and  the  calculated  and  measured  bottom  loss  in  dB  (including  smearing 
effect),  respectively. 

For  the  optimization  of  E  (Ecj.  2)  differential  evolution  (DE)  is  used  ([51).  DE  starts  with 
an  initial  population  of  q  randomly  chosen  parameter  value  combinations  nt.  These  nTs  are 
improved  during  successive  generations  of  constant  size  qy  in  the  sense  that  a  descendant 
replaces  an  m,  becoming  its  successor,  if  it  has  a  lower  energy.  In  [5]  various  ways  to 
generate  descendants,  i.e.,  new  nTs,  are  described*  Here  only  the  following  procedure  is 
adopted.  At  the  start  of  generation  k  the  parameter  vectors  m*,u  ...  mA,<,  are  given  and  for 
each  of  them  a  descendant  is  created,  being  a  potential  successor.  To  create  this  descendant 
dju  of  first  a  partner  p*,,  is  constructed  according  to  =  m* ,/  +  F(m with 
three  different  vectors  m *,,/ ,  m* and  chosen  at  random  from  the  population  and  F  a 

scalar  multiplication  factor  between  0  and  1 .  The  descendant  of  m*,,  results  from  applying 
crossover  to  m*^  and  p*,,  with  crossover  probability  pc*  A  higher  value  of  pc  leads  (on  the 
average)  to  more  dimensions  of  to  be  copied  into  d*.„  while  the  values  for  the  remaining 
dimensions  are  taken  from  m*.,.  Only  if  its  energy  is  lower,  d*,,  replaces  nu„  thereby 
becoming  its  successor.  The  setting  parameters  of  DE  arc  population  size  q ,  multiplication 
factor  F,  crossover  probability  pc  and  the  number  of  generations  Afe.  The  global  search  is 
followed  by  a  local  search  with  the  downhill  simplex  method. 


3.  SUB-BOTTOM  PROFILING  BY  FOURIER  TRANSFORMING  THE 
REFLECTION  COEFFICIENT 

The  reflection  coefficient  is  an  interference  pattern  as  a  function  of  frequency  and  grazing 
angle  caused  by  the  angle  dependent  path  differences  between  the  layer  boundaries.  As 
presented  by  Harrison  ([4])  the  impulse  response  for  each  grazing  angle  is  derived  by  Fourier 
transforming  the  complex  reflection  coefficient.  The  impulse  response  reveals  the  layer 
boundaries  of  the  sea-bottom. 

The  ratio  of  the  downward  and  upward  array  response  gives  the  modulus  square  of  the 
reflection  coefficient  \R\.  However,  the  phase  of  the  reflection  coefficient  can  not  be 
measured.  For  deducing  the  phase,  spectral  factorization  is  applied*  The  impulse  response  r(/) 
now  becomes 


r(0  =  3'jcxp(HH^H)))) 


(3) 


where  3  is  the  Fourier  transform  and  H  is  the  operation  of  adding  the  through  Hilbert 
transformation  generated  imaginary  part  and  original  real  part. 

|  R(o))  |  is  obtained  from  the  data  for  a  limited  range  of  frequencies  only,  with  the  upper 
frequency  (being  1 500  Hz  at  normal  incidence  for  our  data)  and  lower  frequency  limit  caused 
by  the  onset  of  grating  lobes  and  by  the  beam  angle  resolution,  respectively.  Therefore 
|  R(co)  |  is  preprocessed  before  spectral  factorization,  Eq.  (3),  is  used,  extending  the  available 
frequency  range.  |  R  |  is  extrapolated  between  zero  and  the  lower  frequency.  The  extrapolated 
data  between  the  upper  frequency  limit  and  3kHz  will  be  suppressed  by  using  a  windowing 
function.  For  this  extrapolation,  the  modulus  square  reflection  coefficient  is  written  as  a  sum 

of  n  cosines  ([4]),  |fl|  =  £  an  cos(<yr„ ) . 

n=0 

4.  THE  AMBIENT  NOISE  DATA 

For  this  study  we  selected  data  obtained  using  a  drifting  vertical  array  of  hydrophones  in 
the  Mediterranean.  The  array  was  allowed  to  drift  for  around  12  hours  and  covered  in  this 
period  about  5  miles  near  the  Ragusa  Ridge,  south  of  Sicily.  The  experiment  is  described  in 
more  detail  in  [4].  The  array  consisted  of  64  elements  with  an  acoustic  aperture  of  62  m.  Only 
the  middle  32  regularly  spaced  hydrophones  (aperture  16  m,  element  spacing  0.5  m)  were 
used  for  the  processing,  which  consisted  of  beamforming  of  the  noise  data.  The  water  depth 
amounted  to  about  1 30  m.  The  resulting  array  responses  versus  angle  and  frequency  were 
averaged  over  a  series  of  data  spanning  a  few  minutes. 


time  (h) 


Fig.  2:  a)  Inverted  layer  thickness  in  meter  as  a  function  of  time.  The  gray  scale  indicates 
the  corresponding  inverted  sediment  sound  speed  C:  (m/s),  b)  Similar  as  for  a)  but  now  travel 
time. 


5*  RESULTS 


5.1.  Geo- acoustic  parameters  from  inversion  of  ambient  noise  data 

Fig.  2a  depicts  the  layer  depths  as  a  function  of  time  as  obtained  through  inversion  of  the 
ambient  noise  data.  Each  dot  in  the  figure  corresponds  to  one  inversion  result,  he.,  to  one 
measured  reflection  coefficient.  Also  shown  (in  grayscale)  are  the  corresponding  sediment 
sound  speeds.  From  these  results  it  can  be  concluded  that  several  layers  exist  The  sound 
speeds  in  the  area  corresponding  to  the  left  of  the  plot  are  somewhat  higher  than  those  in  the 
area  corresponding  to  the  right  of  the  plot.  In  the  middle  of  the  track  the  layering  shows  a 
very  different  behavior  compared  to  the  layering  at  the  two  ends  of  the  track,  which  is 
probably  caused  by  the  presence  of  a  ridge. 

Fig.  2b  shows  the  same  results  as  Fig.  2a,  but  now  travel  times  (one-way)  (as  derived 
using  the  inverted  C:  values)  are  presented.  Compared  to  Fig.  2a  the  layering  becomes  more 
unique.  It  is  concluded  that  some  of  the  results  of  Fig.  2a  correspond  to  local  optima  that 
differ  in  layer  thickness  and  sediment  sound  speed  but  that  are  equal  with  regards  to  travel 
times.  The  sediment  is  found  to  have  low  attenuation  coefficient,  in  general  less  than  0.2 
d B/X,  with  densities  in  between  1.2  and  1  >8  g/cm\ 


5,2.  Layer  boundaries  from  Fourier  transforming  the  reflection  coefficient 

The  method  of  Fourier  transforming  the  reflection  coefficient  is  applied  to  the  data.  Fig.  3 
shows  the  resulting  profiling  results  constructed  from  the  impulse  responses  between  80°  and 
90°.  Several  layer  boundaries  can  clearly  be  distinguished.  Again  a  change  in  layering  can  be 
seen  at  the  ridge.  The  results  of  Fig.  3  show  good  agreement  with  the  results  presented  in  [4]. 


-0.02  <£ 


’I  0 
§  2 


0  2 
time  (h) 


Fig.  3:  Sub-bottom  profile  derived  with  Fourier  transforming  the  reflection  coefficient.  Layer 
depth  is  in  meters  fusing  a  sediment  sound  speed  of  1500  m/s)  and  array  drifting  time  is  in 
hours , 


53,  Comparison  of  the  tw  o  methods 

Fig.  4  show's  the  comparison  between  the  two  sets  of  estimated  layer  thicknesses.  The  light 
dots  present  layer  thicknesses  as  determined  through  inversion  of  the  ambient  noise  data, 
valid  for  a  sediment  sound  speed  of  1500  m/s.  The  black  dots  indicate  the  first  absolute 


maxima  of  the  inverse  Fourier  transformation  results.  Clearly  the  two  methods  give  layer 
thicknesses  that  are  in  good  agreement. 


*  inverse  Fourier  transformation 
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Fig.  4:  Light  dots:  Inverted  layer  thickness  (for  a  sediment  sound  speed  of  1500  m/s).  Black 
dots:  layer  thickness  obtained  by  the  inverse  Fourier  transformation  method 


6.  CONCLUSIONS 

In  this  paper  two  methods  for  determining  the  seabottom  layering  out  of  ocean  ambient 
noise  have  been  applied  to  a  scries  of  data,  acquired  along  a  5-milc  track.  The  first  method, 
where  the  reflection  coefficient  is  modelled  and  a  global  optimisation  method  applied  to  find 
the  maximum  match  between  modelled  and  measured  data,  docs  not  only  provide  the 
layering,  but  also  allows  for  estimating  other  geo-acoustic  parameters,  such  as  sediment 
sound  speed.  The  computational  load  of  the  method  is  large,  implying  computation  times  of  a 
few  days  on  a  state-of-the-art  PC  The  second  method,  however,  where  the  layering  is 
determined  directly  from  the  data  through  inverse  Fourier  transformation,  gives  seabottom 
layering  only.  The  computational  load  of  this  method  is  small,  resulting  in  computation  times 
of  the  order  of  minutes. 

By  comparison  of  the  two  sets  of  layer  estimates  it  has  been  demonstrated  that  both 
methods  give  results  that  are  in  good  agreement  and  apparently  perform  equally  well  for  sub¬ 
bottom  profiling.  If  one  is  interested  in  sediment  classification,  the  inverse  Fourier  transform 
method  could  be  used  as  a  pre-processing  step  for  the  inversion,  providing  the  number  of 
layers  and  their  order  and  an  estimate  for  their  thicknesses. 
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